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Many complex systems exhibit hydrodynamic (or macroscopic) behavior at large scales character- 
ized by few variables such as the particle number density, temperature and pressure obeying a set of 
hydrodynamic (or macroscopic) equations. Does the hydrodynamic description exist also for waves 
in complex open media? This is a long-standing fundamental problem in studies on wave localiza- 
tion. Practically, if it does exist, owing to its simplicity macroscopic equations can be mastered far 
more easily than sophisticated microscopic theories of wave localization especially for experimental- 
ists. The purposes of the present paper are two-fold. On the one hand, it is devoted to a review of 
substantial recent progress in this subject. We show that in random open media the wave energy 
density obeys a highly unconventional macroscopic diffusion equation at scales much larger than 
the elastic mean free path. The diffusion coefficient is inhomogeneous in space; most strikingly, as 
a function of the distance to the interface, it displays novel single parameter scaling which captures 
the impact of rare high-transmission states that dominate long-time transport of localized waves. 
We review aspects of this novel macroscopic diffusive phenomenon. On the other hand, it is devoted 
to a review of the supersymmetric field theory of light localization in open media. In particular, 
we review its application in establishing a microscopic theory of the aforementioned unconventional 
diffusive phenomenon. 

PACS numbers: 42.25.Dd, 71.23.An 

I. INTRODUCTION AND MOTIVATION 

Anderson localization is one of the most profound concepts in modern condensed matter physics (see Refs. [T] and 
[2] for recent reviews). While this phenomenon was originally predicted for electron systems [3j, direct observations of 
electron wave localization are notably difficult because of electron-electron interactions. In the eighties, the universality 
of Anderson localization as a wave phenomenon was appreciated by many researchers [IHSj- The study of classical 
wave localization is now a flourishing field [5Hll). In recent years the unprecedented level reached in manipulating 
dielectric materials [T2H22] and elastic media [231 El] has led to substantial experimental progress in studies of 
Anderson localization in various classical wave systems. Together with the realization of Anderson localization in 
ultracold atomic gases |251 126) . these experiments have raised many fundamental issues, and are triggering a renewal 
of localization studies. In addition, the experimental observation of Cao and co-workers [17] has triggered considerable 
investigations of using Anderson localization of light for random lasing (see Ref. [28] for a review) . 

Since photons do not mutually interact, it was anticipated [S1IS1[5S] that classical electromagnetic waves may serve 
as an ideal system for experimental studies of Anderson localization. Ground-breaking experimental achievement was 
made in 2000 when unambiguous evidence of microwave localization in quasi one-dimensional (QID) samples was 
observed I15J. Since then substantial experimental progress in electromagnetic wave localization have been achieved. 
Among them are dynamics of localized microwave radiation in QID samples |21U22| . time-resolved transmission of light 
through slab media [T^ HD] , measurements of the spatial distribution of the localized modes [5D] , and the observation 
of two-dimensional Anderson localization in photonic lattices [T7]. It has also been within the reach of microwave 
experiments on QID localized samples the statistics of quasi-normal modes |31j and transmission eigenvalues |32j . 
In particular, the crystallization of transmission eigenvalue distribution [33-38] , has been observed recently. While 
these experimental studies have not yet been carried out in other wave systems, they further show that classical 
electromagnetic waves have great advantages in experimental studies of localization. 

Importantly, experimental setups for probing localization of electromagnetic waves [T5| 1191 [21] [22l l39l - l42| and other 
classical waves [21j typically are very different from those for quantum matter waves [25l, . In the former, wave 
energies leak out of random media through interfaces and measurements are performed outside media. Specifically, 
in transmission experiments, waves are launched into the system on one side and detected on the other (Fig. [T]); in 
coherent backscattering experiments (see Ref. j43j for a review), the light beam is launched into and exit from the 
medium at the same air-medium interface, and the angular variation of the reflected intensity is measured (Fig. [2|. 
Theoretically, one treats the former system as a finite-sized open medium and the latter one (often) as a semi-infinite 
open medium. Therefore, these classical wave experiments address a fundamental issue - the localization property of 
open systems - which, as we will see throughout this review, conceptually differs from that of bulk (infinite) systems. 

In fact, open systems have been at the core of localization studies for several decades. In the mid-seventies, Thouless 
and co-workers noted that the unique parameter governing the evolution of eigenstates when the system's size is scaled 
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FIG. 1: Example of light transport through open media. Upper panel: microwave radiation is launched from a horn placed 
before the sample, a copper tube consisting of randomly positioned alumina spheres. Wave intensity on the output plane is 
measured by the detector placed in front of the sample. Lower panel: intensity spectra at different positions ri.2 of the output 
plane. (Adapted by permission from Macmillan Publishers Ltd: Nature 471, 345-348, © 2011) 



is the dimensionless conductance (which is nowadays called the Thouless conductance), a main transport characteristic 
of open (electronic) systems [44 . This eventually led to the advent of the milestone discovery ~ the single parameter 
scaling theory of Anderson localization - in 1979 [IS]. Shortly later, Anderson and co-workers further pointed out 
[3S] that since the resistance displays a broad distribution in the (one-dimensional) strongly localized regime, it is 
not sufficient to study the scaling behavior of the average conductance. Instead, one must study how the entire 
conductance distribution evolves as the sample size increases. In fact, long before the Anderson localization theory 
was developed, such a distribution was discovered by Gertsenshtein and Vasil'ev in a study of the exponential decay 
of radio waves transporting through waveguides with random inhomogeneities. Importantly, the large-conductance 
tail of the distribution represents rare localized states peaked near the sample center with exponentially long lifetime 
P8 51 . These states have been found experimentally to be responsible for long-time transport of localized waves 
through open systems |21| . They have high transmission values that may be close to unity, in sharp contrast to typical 
localized states with exponentially small transmission. Such peculiar features are intrinsic to open medium. The high 
transmission states even promise to have practical applications: they mimic a 'resonator' with high-quality factors 
(due to the long lifetime) in optics and thus can be used to fabricate a random laser |52j . In combination with optical 
nonlinearity, it can also be used to realize optical bistability |53| . 

Nowadays experimental and theoretical results on global transport properties (in the sense of that they provide no 
information on wave propagation inside the medium) - characterized by conductance, reflection, transmission, etc. 
- of both classical and de Broglie waves in one-dimensional open media have been well documented (for a review, 
see, e.g., Ref. [H]). An essential difference between finite-sized open system and infinite (closed) system was pointed 
out by Pnini and Shapiro [56] . That is, the wave field is a sum of traveling waves for the former and of standing 
waves for the latter. Let us mention a few more recent results in order for readers to better appreciate the rich 
physics arising from the interplay between localization and openness of the medium. In Ref. |18j . it was found that, 
surprisingly, even low absorptions essentially improve the conditions for the detection of disordered-induced resonances 
in reflection as compared with the absorptionless case. In Ref. [57 , Fyodorov considered reflection of waves injected 
into a disordered medium via a single channel waveguide, and discovered the spontaneous breakdown of S'-matrix 
unitarity. Interestingly, this may serve as a new signature of Anderson transition in high dimension. Much less 
is known regarding how wave propagates inside the medium. Nevertheless, they may provide a key to many new 
experimental results such as the spatial distribution of localized modes [30] and the interplay between absorption 
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(gain) and resonant states [15] . 

In principle, one may describe the wave field in terms of the superposition of quasi-normal modes \21\ I58j 
namely the eigenmodes of the Maxwell equation in the presence of open boundary. This exact approach carries full 
information on wave propagation but in general must be implemented by experiments or numerical simulations. As we 
are interested in physics at scales much larger than the elastic mean free path, alternative yet simpler approach might 
exist. The situation might be similar to that of many complex systems like fluids. There, many degrees of freedom 
causing the complexity of dynamics notwithstanding, at macroscopic scales (much larger than the mean free path) the 
system's physics is well described by few variables such as the particle number density, temperature and pressure and a 
set of hydrodynamic (or macroscopic) equations. In fact, for open diffusive media, such a macroscopic equation is well 
known which is the normal diffusion equation of wave energy density reflecting Brownian motion of classical photons. 
For open localized media, the normal diffusion equation breaks down. The question of fundamental interests and 
practical importance ~ addressed from experimental viewpoints by many researchers for decades [21, |i23, 43, 64-65] 
- thereby arises: is the macroscopic description valid for open localized media? Since the eighties j61ti63j there have 
been substantial efforts in searching for a generalized macroscopic diffusive model capable of describing propagation 
of localized waves in open media. In the past decade, research activity in this subject has intensified pTl [M[ IMHTOj . 
In essence, one looks for certain generalization of Pick's law with the detailed knowledge of multiple wave scattering 
entering into the generalized diffusion coefficient. Such a diffusive model has the advantage of technical simplicity over 
many sophisticated first-principles approaches, and may provide a simple principle for guiding experimental studies 
(see Ref. [43] for a review). 

Important progress was achieved by van Tiggelen and co-workers in 2000 [53]. These authors noticed that in 
open media weak localization effects of waves are inhomogeneous in space. Therefore, they hypothesized the position 
dependence of the one- loop weak localization correction to the diffusion coefficient. By further introducing a key 
assumption - the validity of the one-loop self-consistency, they obtained a phenomenological, nonlinear macroscopic 
diffusion model describing static transport of localized waves in open media. Later, the proposed self-consistent local 
diffusion (SCLD) model was extended to the dynamic case [651I66]. In the past decade the SCLD model has been 
used to guide considerable research activities in classical wave localization. The subjects include: static reflection 
and transmission of light in systems near the mobility edge |64| . dynamic reflection and transmission of QID weakly 
localized waves [65| . dynamics of Anderson localization in three-dimensional media |66j . transmission of localized 
acoustic waves through strongly scattering plates ^24 , and transmission and energy storage in random media with 
gain [71]. 

This prevailing model was first examined in Ref. [21| in experiments on pulsed microwave transmission in QID 
localized samples. The SCLD model failed to describe the experimental results. In fact, in these systems, the long- 
time transport of localized waves is determined by rare disorder-induced resonant transmissions. In addition, the 
SCLD model fails to describe results obtained from numerical simulations on static (long-time limit) wave transport 
through one-dimensional localized samples [BHl [7D] . Since the macroscopic diffusive model describes local (in space) 
behavior of wave propagation (in the average sense), it is very 'unlikely' that it could capture simultaneously two 
prominent features of resonant transmissions. That is, they are rare events and objects highly non-local in space. 
In fact, it had been commonly suspected that eventually the validity of the macroscopic diffusion concept in open 
systems would be washed out by these rare events. On the other hand, the predictions of the SCLD model agree with 
the experimental results in higher dimensions [24] . In view of these contradictory observations one inevitably has to 
re-examine the fundamental issue: whether and how do localized waves in open media exhibit macroscopic diffusion? 

A preliminary but very surprising answer to this question was provided in a first-principles study of static transport 
of localized waves in one-dimensional open media [69] . It turns out that in these systems, the macroscopic diffu- 
sion (or equivalently Pick's law) concept never breaks down at scales much larger than the elastic mean free path. 
Rather, in depending on the distance to the interface, the diffusion coefficient exhibits novel scaling. It is such scaling 
that unifies the objects of resonant transmission and macroscopic diffusion. The novel scaling behavior is missed 
by the phenomenological SCLD model. In this review, we will discuss the existence of such highly unconventional 
macroscopic diffusion in more general situations (e.g., dynamic transport and in high dimensions). It is believed that 
waves propagating through open media may exhibit even richer macroscopic behavior. As such, the unconventional 
macroscopic diffusion of waves in open media may potentially open a new direction for the study of Anderson local- 
ization. One of the main purposes of this paper is to review various aspects of this diffusion phenomenon such as its 
microscopic mechanism, microscopic theory and perspectives. 

Technically, treatments of localization in open media differ dramatically from those of infinite media. The study 
of open systems has proved to be an extremely difficult task especially for high dimensions and in the dynamic 
case. Recall that for infinite electron systems, the diagrammatic technique [75H75] has had considerable successes in 
(approximate) studies of a variety of localization phenomena, e.g., one-loop weak localization, strong localization, and 
criticality of Anderson transition (see Refs. [TIKT^ for reviews). In view of these successes, it is natural to extend this 
technique to classical wave systems, and this has been done by many authors [77H75] . In the past, the non-perturbative 
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FIG. 2: Example of light transport through open media. Experimental setup (upper panel) of coherent backscattering (lower 
panel): a beam of light is launched into and backscattered by a semi-infinite random medium, (upper panel: from Ref. [109j 
with reproduction permission from G. Maret © The American Physical Society) 



diagrammatic theory of Vollhardt and Wolfle (VW) [751 [73] , having the advantage of technical simplicity, has become 
a popular approach in studies of classical wave localization. Furthermore, it has been noticed (MHSSl UHl [HHI that it is 
possible to generalize the VW theory appropriately to describe transport of strongly localized waves in open media. 
It should be emphasized that the rational of the original VW theory [72j [73] is built upon the exact Ward identity 
[7?l 1751 [5T| and the sophisticate summation over the dominant infrared divergent diagrams (see Refs. [TH [55] for a 
technical review). This rational is at the root of the one-loop self-consistency and to the best of my knowledge, has 
not yet been established for open media. These considerations necessitate the invention of an exact theory beyond 
simple one-loop perturbation. 

There have been few attempts |3H71 [53] of generalizing the replica field theory [5^H55] and the Keldysh field theory 
[57] to classical wave systems. However, these studies [3H3 dH] focus on infinite random media. For open systems 
some exact solutions for QID strong localization have been found by using the replica field theory [38]. As for the 
replica field theory, in general, its applicability in the strongly localized regime needs further investigations due to 
the well-known problem of analytic continuation [55] [5^ . As for the Keldysh field theory, it does not encounter such 
a difficulty. However, to the best of my knowledge, it is not clear how to use this theory to obtain concrete results 
for strong localization. Finally, the Dorokhov-Mello-Pereyra-Kumar (DMPK) equation [33 1 [34 [ [90 ] is an exact theory 
and a very powerful approach to QID strong localization in open media. The special case of the DMPK equation 
in one dimension was discovered as early as in 1959 |47| . But, it mainly provides information on global transport 
properties such as transmission and its fluctuations; it provides no information on local wave transport properties: it 
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could not describe how waves propagate from one point to the other inside media. For example, by using this theory 
one could not judge the possibility of generalizing Fick's law to localized open media. In addition, it is well-known 
that the DMPK equation is valid only in one dimension. 

The supersymmetric approach escapes all these key difficulties. In Ref. [91 , the supersymmetric quantum me- 
chanics was used to study reflection of classical waves by one-dimensional semi-infinite disordered medium. The 
supersymmetric field theory invented by Efetov [88l |92l |93j was used by Mirlin and co-workers to study the deviation 
from the Rayleigh distribution of classical light intensity in diffusive QID samples (see Ref. [Ml for a review). In 
Ref. [67j , the supersymmetric field theory was generalized to high-dimensional open media with internal reflection, 
and was employed to investigate dynamic transport of classical waves in high-dimensional open media. In particular, 
the unconventional macroscopic diffusion was discovered in Ref. |69j from this first-principles theory |67j . In fact, the 
machinery of supersymmetric field theory has become a standard approach in studies of disordered electronic systems. 
Nonetheless, researchers working on classical wave localization are less familiar with this technique. On the other 
hand, there are a number of recently observed phenomena regarding classical wave transport through open media 
which can be thoroughly studied by using this technique. These include effects of internal reflection on transmission 
fluctuations '32] and dynamics of localized waves [21' . Since there are many prominent differences between classical 
and electronic waves (e.g., the condition of strong localization [Hd]), and a review of the application of the super- 
symmetric field theory in light localization has been absent, it is necessary to introduce here this technique - in the 
context of transport of classical electromagnetic waves through open random media - at a pedagogical level, which 
serves another main topic of this review. 

The review is written in a self-contained manner, with the hope that readers wishing to master the supersymmetric 
technique and then apply it to classical wave systems could follow most technical details without resorting to further 
technical papers. Because of this, it could not be a complete introduction of this advanced theory. In fact, there have 
been a comprehensive book [88' and several excellent reviews pi l54[[?5ll93) covering various aspects of supersymmetric 
field theory. In particular, we will not introduce non-perturbative treatments [371 f53ti55] , because these are technically 
highly demanding. The present review, implemented with substantial technical details, may be considered as a 'first 
course in practicing the supersymmetric field theory'. It aims at helping the researchers working in the field of 
classical wave localization to become familiar with formulating large-scale wave dynamics in terms of supersymmetric 
functional integral formalism and further to master its perturbative treatments. With these preparations, readers 
may be ready to manipulate its non-perturbative treatments after further resorting to original technical papers and 
reviews. We emphasize that although the present review is written in the context of classical waves, unconventional 
macroscopic diffusion to be reviewed below is a universal wave phenomenon in open systems. In particular, it also 
exists for de Broglie waves. 

For simplicity, we shall focus on scalar waves throughout this review. We shall not distinguish the concepts of 
(classical) scalar wave, electromagnetic wave, and light. The remainder of this review is organized as follows. In 
Sec. [Ill we will review a number of research activities in studies of transport of localized waves through open media. 
These studies may be roughly classified into two categories: one is based on the macroscopic diffusion picture and the 
other on the mode picture. In Sec. |III| we will discuss in a heuristic manner the 'origin' of supersymmetry in studies of 
disordered systems. In Sec. |IV| we will proceed to introduce the supersymmetric field theory suitable for calculating 
various physical observables in random open media. Then, in Sec. [Vj we will apply this theory to study a special 
quantity, the spatial correlation of wave intensity. In particular, we will develop a general theory of unconventional 
macroscopic diffusion for the general case of a random slab. In Sec. |VI[ we will study in details the special case of 
one-dimensional unconventional macroscopic diffusion. In particular, we will present the explicit results for the static 



local diffusion coefficient and provide corroborating numerical evidence. Finally, we conclude in Sec. VII In order to 
make this review self-contained, we include a number of technical details in Appendices [A|e1 



II. MACROSCOPIC DIFFUSION VERSUS MODE PICTURE OF WAVE PROPAGATION 

In infinite three (or higher)-dimensional systems, the envelope of eigenfunctions decays exponentially in space for the 
distance from the center far exceeding the localization length, ^, when disorder is strong. Below or in two dimensions, 
even weak disorders lead to localization 15^, 'M', 1 55141011 [lll| . The (asymptotic) exponential decay of eigenfunctions 
in space is a key characteristic of localized waves in the mode picture. Can localized waves be understood in terms 
of the macroscopic diffusion picture? Recall that in the absence of wave interference effects, light exhibits pure 
'particle (classical photons)' behavior, which in random environments is the well-known Brownian motion (We shall 
not consider here the case where the motion of photons is non-Brownian.) or normal diffusion of photons. When wave 
interference is switched on, localization effects set in. It has been well established that the latter does not render the 
macroscopic diffusion concept invalid [7^H751 [551 IMl I100H102| . Rather, they (strongly) renormalize the Boltzmann 
diffusion constant. Specifically, given the point like source Jt(r') located at r', the disorder averaged intensity profile 
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at time t is given by [51] /(r,t) = / |^ |^e*"*>'(r, r'; w) J^^l (r') J*_ji (r'), where Jui{r') is the spectral decomposition 
of the source. Here, we have introduced the spatial correlation function defined as 

in the frequency domain. Here, G^''^^, are the retarded (advanced) Green functions whose explicit definitions 



will be given in Sec. HI A Throughout this review we use (•) to denote the disorder average. Most importantly, the 



correlation function satisfies [731 EH] 

(-iw - D{u)\'^)y{r, r'; w) = S{r - r'). (2) 

It is important that here, the diffusion coefficient, D{lo), is frequency-dependent. The latter accounts for the retarded 
nature of the response to the spatial inhomogeneity of the wave energy density. 

How are the mode and macroscopic diffusion pictures modified in open media? For the mode picture, because of 
the energy leakage through the system's boundary, the Hermitian property is destroyed. Consequently, one describes 
transport in terms of quasi- normal modes [1TJ[3I1[SH] each of which has a finite lifetime. Loosely speaking, wave energy 
is pumped into these modes, stored there and emitted in the course of time. For the macroscopic diffusion picture, 
substantial conceptual issues arise. In this section, we will review great efforts made in extending the macroscopic 
diffusion concept to localized wave transport through open media. 



A. Electromagnetic wave propagation: macroscopic diffusion picture 

Macroscopic diffusive approach to light transport may be dated back to the first half of last century, when the 
Boltzmann-type kinetic equation was used to study radiation transfer in certain astrophysical processes [103] . In fact, 
it is a canonical procedure of deriving a normal diffusion equation - valid on the scale much larger than the transport 
mean free path, I - from the Boltzmann kinetic equation. In the normal diffusion equation, the diffusion constant 
D{uj) = cl/d = Dq, where d is the dimension and c the wave velocity. Notice that the latter may be reduced by the 
scattering resonance |104j which we shall not further discuss, and from now on we set c to unity. 

The Boltzmann transport theory essentially treats light scattering off dielectric fiuctuations as Brownian motion of 
classical particles, the 'photons', and ignore the wave nature of light. However, the latter gives rise to interference 
effects which have far-reaching consequences. A canonical example is the coherent backscattering of light from semi- 
infinite diffusive media (Fig. [2] lower panel). There, two photon paths following Brownian motion may counterprop- 
agate and thus constructively interfere with each other jl05| I106j , leading to enhanced backscattering [43l llOTHllOj . 
Effects of wave interference are even more pronounced in low dimensions (d < 2) or in three-dimensional strongly 
scattering random media, where light localization eventually occurs [5HSj- To better understand the macroscopic 
diffusion picture of wave transport in open localized media we begin with a brief summary of the (conventional) 
macroscopics of wave propagation in infinite media. 



1. Infinite random media 

In this case, the diffusion of wave energy follows Eq. The leading order wave interference correction to Da - 
the well-known weak localization correction due to constructive interference between two counterpropagating paths - 
was first found (for electronic systems) using the diagrammatic method [721 1731 176] , which is 

SD^>{uj)^ / T^r-u -,012 ' (3) 

TTi/ J {2tt)"- —lUJ + Doii 

with v is the local density of states of particles. The result is perturbative and valid only if SD^^^Q)/ Dq ^1. In 
spite of this limitation, Eq. ([s]) has many important implications: the lower critical dimension of Anderson transition 
is two, because it suffers infrared divergence for dimension d < 2 [JS]. That is, the system is always localized in low 
dimensions {d < 2) while exhibits a metal- insulator transition in higher dimensions {d > 2). 

(i) In dimension d < 2, ^Z?^^)(<D) is much smaller than Dq for large frequencies, ui 3> Dq/C^- This implies that at 
short times (^ ^'^/Dq), wave transport largely follows normal diffusion with a small suppression due to interference. 
For low frequencies, ui ^ Dq/^^ which corresponds to long times (^ ^'^/Dq), one finds using non-perturbative 
methods [72. ,73. .94. .lOUl fTUTl llll| that the dynamic diffusion coefficient, D[Cj), crosses over to ~ —iCoS,'^. Notice that 
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is a hallmark of strong localization. Namely, wave energy diffusion stops, (ii) In dimension d > 2, result 
([3| does not diverge in the infrared limit {uj — 0) but suffers an ultraviolet divergence. To cure the latter, one needs 
to introduce the ultraviolet cutoff ~ l^^ for the integral over k. This introduces a dimensionless coupling constant 
vDqI'^~^ ^ Dq/SD^^\lli = 0). For large coupling constant, the interference correction ^ is negligible; for a coupling 
constant order of unity, the interference correction and Dq are comparable. The latter is the onset of the Anderson 
transition, and go = i'DqI'^^^ — 0{1) is equivalent to the well-known loffe-Regel criterion |112| . Above the critical 
point, the system is in the metallic phase and D(u)) approaches some non-zero constant in the limit J) — > 0; below the 
critical point, the system is in the localized phase and D{oj) vanishes again as ~ — ztD^^. (iii) At the critical point (in 

dimension d > 2), the low- frequency behavior of the diffusion coefficient is given 

These results have also been predicted for completely different systems - the quantum kicked rotor jll41 1115j - by 
using first-principles theories |116l 1117] , and have been observed experimentally |118H120| . 

Let us make two remarks for the macroscopic diffusion equation First, the dynamic diffusion coefficient D{Cj) 
is homogeneous in space. Secondly, this shows that wave interference does not destroy Fick's law: the energy flux is 
proportional to the (negative) energy density gradient, given by — Z?((I')V/(r; w) in the frequency domain. According 
to this law, the response to inhomogeneous wave energy density is local in space and non-local in time (retarded 
effect). 



2. Open random medta 

The study of light propagation in open random media has a long history [103j , and since the eighties researches on 
this subject have been strongly motivated by the search for light localization. As mentioned above, in the diffusive 
regime, wave interference corrections to the Boltzmann diffusion constant are negligible and photons follow normal 
diffusion or Brownian motion. The effect of the air-medium interface is to introduce certain boundary conditions 
that implement the normal diffusion equation (see Ref. |104j for a review). The boundary condition is essential to 
calculations of coherent backscattering lineshape (for diffusive samples) |1071I110| . In this case, the boundary condition 
effectively enlarges the diffusive sample, and the effective interface is located outside the medium at a distance oi C, - 
the so-called extrapolation length - to the genuine interface. The relation between the extrapolation length and the 
internal reflection coefflcient was studied by various authors |I21HI23] . 

Can wave propagation in open localized media he described by certain macroscopic diffusion equation? This funda- 
mental problem has attracted the attention of many researchers. The flrst attacks were undertaken more than two 
decades ago (e.g., Refs. [6Tti63] ). motivated by experiments on coherent backscattering of light from strongly disordered 
media. The earlier attempts resort to phenomenological generalization of the scale-dependent diffusion coefficient de- 
veloped for infinite media to open media |I24Uf25] . While these theories had been debated, an important observation 
was made in Ref. [64 . That is, because photons near the air-medium interface easily escape from the medium, the 
returning probability density near the interface must be smaller than that deep in the medium. Consequently, wave 
interference (or localization) effects must be inhomogeneous in space. The strength of interference effects increases 
as waves penetrate into the medium. Indeed, both the supersymmetric field theory 'BT and the diagrammatic theory 
[68j developed very recently for open media justify this important observation. Specifically, it is shown that near the 
interface the diffusion coefficient may be largely unrenormalized even though strong localization develops deep in the 



medium (cf. Sec. VI I. Interestingly, a similar conclusion was reached in the study of superconductor- normal metal 
hybrid structures even earlier |126j . 

To take the inhomogeneity of wave interference effects into account, the phenomenological SCLD model was intro- 
duced in a series of papers |64H66j which include two key assumptions. The first is to hypothesize a local (position- 
dependent) diffusion coefficient, D{y;Cj), that locally generates a macroscopic energy flux via Fick's law. Specifically, 
instead of Eq. ([2]), the correlation function y{r,r';Cj) satisfies 

[-iCo - V • D{v, Cj)V] y{v, r'; w) = 5{r - v'). (4) 

The crucial difference from Eq. ^ is the position-dependence of the diffusion coefficient. The second is the one- 
loop consistency which is a phenomenological generalization of the VW theory [751 [73]. That is, the local diffusion 
coefficient is given by 



D{v-w) Do 



(5) 



It was argued that the SCLD model is valid for both weakly and strongly disordered systems. As mentioned 

in the introductory part, the SCLD model has guided a number of research activities on localization in open media. 
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FIG. 3: A Gaussian pulse is launched into QID samples of different sample lengths (61 and 90 cm), and the averaged intensity 
{I{t)} and decay rate {A{t)} on the output plane are measured. The experimental results (solid, in black) are compared with 
the predictions of the SOLD (namely 'SCLT' in panel (a)) model (dashed, in red) and of classical diffusion theory (dotted, in 
blue). Although the prediction of the SCLD model and the measurements agree well at short times ((a)), dramatic deviations 
show up at long times ((b)-(d)). The curves are normalized to the peak value in (a) and (b). (from Ref. [2T] with reproduction 
permission from Z. Q. Zhang © The American Physical Society) 



In particular, the prediction of the SCLD model and measurements in experiments on localized elastic waves in three 
dimensions agree well [24_. 

In (quasi) one-dimension strong disagreement between the predictions of the SCLD model and measurements in 
experiments and numerical simulations has been seen [211 169| . In the dynamic case, cD 7^ 0, both experiments and 
numerical simulations on the dynamics of microwave pulse propagation through QID samples have been carried out 
[3T], and the results of the time-resolved transmission are compared with the prediction of Eqs. ^ and ([s]). As shown 
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in Fig. [3] although the SCLD model can account very well for the time-resolved transmission at intermediate times, 
it fails completely at longer times. In the steady state {Co = 0), numerical simulations of the wave intensity profile 
across the sample have been carried out and the local diffusion coefficient was computed [53]. The numerical results 
are compared with the predictions of the SCLD model, and dramatic deviations are found (see Fig.|4|. 

It has been further shown |21l 169] that in (quasi) one-dimensional systems transport of localized waves is dominated 



by rare disorder-induced resonant transmissions which lead to far-reaching consequences (see also Sec. VI I. In fact, 
by using the first-principles microscopic theory it has been shown [69J that for x deep inside the samples, i.e., 
^ <C min(x, L — x), 

x(I,-^) 

D{x)^e (6) 

This result has been fully confirmed by numerical simulations (see Fig.|4]). Notice that this expression is symmetric with 
respect to the sample mid-point, i.e., D(x) = D{L — x). Most importantly, near the sample center, the enhancement 

(min(j;, 1,-x))^ I I 

from the prediction of the SCLD model is exponentially large ^ e ^« . As to be shown in Sec. VI this dramatic 
enhancement finds its origin at the novel scaling displayed by D{x) [53]. That is, D{x) depends on x via the scaling 
factor x{L — x)/{L^). 




10 
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FIG. 4: Results of D{x)/Do obtained from numerical simulations (squares and circles), the analytic prediction Q (solid lines) 
and the SCLD model (dashed lines) are compared. Numerical simulations are performed for two (angular) wave frequencies, 
tj = 1.65c/a (square) and oj = 0.72c/a (circle), and for five different sample lengths, L/^ — 2.5, 5, 10, 15 and 20. (from Ref. |69] 
with reproduction permission from C. S. Tian, S. K. Cheung, and Z. Q. Zhang © The American Physical Society) 



3. Is the concept of local dijfusion universal? 

In some optical systems (e.g., Faraday- active medium [1071 1129] ) the one-loop weak localization may be strongly 
suppressed and eventually time-reversal symmetry may be broken. In these cases, the SCLD model is no longer 
applicable since it crucially relies on the one-loop self-consistency or the time-reversal symmetry. On the other hand, 
there have been rigorous studies showing that systems with/without the time-reversal symmetry (more precisely, 
corresponding to the Gaussian orthogonal/unitary ensemble (GOE/GUE) in the random matrix theory [SU]) have 
largely the same strong localization behavior, and the only difference is the numerical factor of the localization length 
[94L 1951 197L I98j . An important question therefore arises: Is local diffusion an intrinsic macroscopic phenomenon of 
(localized) open media? This question was first studied in Ref. |69j . By using the first-principles theory to be reviewed 
below the authors showed that Eq. ^ is valid for both GOE and CUE systems, and the symmetry only affects the 
numerical coefficient of the localization length, ^. 



B. Electromagnetic wave propagation: mode picture 



Alternatively, transport of waves through open media may be understood in terms of quasi-normal modes. Each 
quasi-normal mode (labeled by n) is characterized by the (complex) wave amplitude, (^„(r), the central frequency a;„, 
and linewidth r„ > 0. The last represents the rate of wave energy leakage through the interfaces. Upon pumping 
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wave energies into the system, a number of quasi-normal modes are excited, and the wave energies are stored in 
these modes. In terms of the quasi-normal mode picture, the essential differences between diffusive and localized 
samples are as follows. For the former system the quasi-normal modes are extended in space, and wave energies are 
readily transported through the sample. Correspondingly, the lifetimes of quasi-normal modes are short, resulting in 
a characteristic linewidth greatly exceeding the mean spacing between neighboring eigenfrequencies. For the latter 
system, there are quasi-normal modes j48H50j whose coupling to outside environments is exponentially small with 
linewidths much smaller than the mean eigenfrequency spacing. At long times, wave energies are stored mainly in 
these long-lived modes which therefore play decisive roles in transport of localized waves. 



1. Dynamic single parameter scaling model 



Dynamics of localized waves undergoes an essential change in the modal distribution of energy with time: the 
transmission is mainly due to short-lived overlapping modes at early times and by long-lived localized states at longer 
times. The long-time dynamics of localized waves is dominated by the spectrally isolated and long-lived localized 
modes, and is well captured by a phenomenological model, the so-called dynamic single parameter scaling (DSPS) 
model |2TJ im , developed for (quasi) one-dimensional systems. 

Consider a one-dimensional sample of length L 3> ^ with transparent interfaces. For a resonantly excited localized 
mode peaked at a distance x from, say, the left interface (x — 0), the steady-state intensity at the right interface 
{x = L) relative to the incident wave is [151 [50] . 



-l(L-2x) 



< X < L/2. 



(7) 



Notice that is symmetric with respect to the sample center, i.e., = T^-x- According to the one-parameter scaling 
hypothesis [15], the Lyapunov exponent (namely the inverse localization length) 7 follows the Gaussian distribution. 



^(7) 



exp 



1\2 



(7-r^) 



(8) 



For (spectrally isolated and long-lived) localized states, the decay rate is the ratio of the sum of the outgoing fluxes 
at two interfaces to the integrated wave energy inside the sample. 



r(7,a;) - 7 



1 



(9) 



Here, an overall numerical factor has been ignored since it is irrelevant for further discussions. Assuming that the 
modes are uniformly distributed inside the sample, we find the transmitted intensity to be 



I{t) = I{x = L,t) 



Ji/L Jo 



(10) 



The theoretical prediction of the DSPS model (Eqs. ([s]) and (9])) for the decay rate —d\nl{t)/dt is in good agreement 
with both the experimental measurements and the numerical simulations [21]. 



2. Complete modal analysis 

Complete modal analysis has been further performed experimentally by Wang and Genack [31j . The underlying 
general principle is as follows. (In fact, a similar principle has been adopted in studies of transport through quantum 
chaotic systems, see Ref. [116] for a review.) For a linear medium, the wave field, denoted as E{r, t), is the superposition 
of quasi-normal modes, 

i?(r, t)=Y. C"^" W e-^(""-^r")*, (11) 

n 

where the coefficients, Cn's, are fixed by the initial condition. The wave intensity profile is given by 

\E{v,t)\' = ^ |c„^„(r)|2e-r"* + ^ c„c;,(p„(r)v,„, (r)*e-*(""-""')*e-(r.+r„0t/2. (12) 



11 



The first term is an incoherent contribution reflecting wave transport through individual modes. In contrast, the 
second term introduces interference between different modes. Since the system is finite, the central frequencies are 



discrete with a characteristic spacing Auj. At long times, t l/Aw, the second term becomes negligible and Eq. (12) 
is simplified to 

\E{r,t)\' ^ |c„^„(r)p e'^-'. (13) 

n 

This shows that transport at long times is dominated by the long-lived modes with small r„. 

It is not an analytic tractable task to obtain complete knowledge on the spectrum and the wave amplitude. Wang 
and Genack realized |31| that quantitative analysis may be substantially simplified in combination with experimental 
measurements. Experimentally, a pulse is incident from one interface and propagates to the other. The field pattern at 



the output surface is recorded, which depends only on the coordinates in the transverse plane. According to Eq. ( 11 ), it 
can be decomposed in terms of the so-called volume field speckle pattern, i.e., c„(/5„(r) with the longitudinal coordinate 
fixed at L. Then, both the spectrum and the corresponding mode speckle pattern can be determined experimentally, 
see Fig. [ij lower panel for typical intensity spectra which is the squared modular of mode speckle pattern. They in 
principle afford a full account of dynamic and static transmission. Experiments confirmed a broad range of decay 
rates and that at long times, the transmission is indeed dominated by incoherent contributions [31] . 



III. 'ORIGIN' OF SUPERSYMMETRY 



How does the supersymmetry enter into the theory of Anderson localization? In this section we will discuss in a 
heuristic manner how the supersymmetric trick is prompted in studies of disordered systems particularly in Anderson 
localization systems [127( 1128] . Furthermore, we will present some crude technical hints showing how the use of this 
trick eventually leads us to a supcrmatrix field theory. The introduction of the supermatrix field lies at the core 
of nonlinear supermatrix a model theory of localization, as will become clearer in the remainder of this review. To 
understand the content of this section the readers need to have some basic knowledge on the Grassmann algebra and 
supermathematics, a preliminary introduction to which is given in Appendix [A] 



A. The supersymmetric trick 



Consider a random medium embedded in the air background. Microscopically, the wave field, E, is described by 
the Hclmholtz equation [51 [71 1130] 



{V2+c^2(l + e(r))}£;(r) =0, 



(14) 



where lu is the (circular) frequency. This equation also describes the propagation of elastic waves [11 [5l [9l [79] , and is a 
good approximation [TDl] to electromagnetic waves provided that the vector character (which is essential, for example, 
to multiple scattering in a Faraday- active medium in the presence of a magnetic field [13111131? ) is unimportant. e(r) is 



the fluctuating dielectric field (with zero mean). Interestingly, classical scalar wave equation (14 1 bears a firm analogy 
to the Schrodinger equation: H = — — a;^e(r) plays the role of the 'Hamiltonian' and the 'particle energy'. 
Notice that here the 'potential', a;^e(r), is 'energy'-dependent. (The interesting property of the energy-dependence of 
the potential immediately leads to important consequences of light localization which are profoundly different from 
electron localization. That is, in the zero frequency limit, the system has an extended state where the wave field is 
uniform in space. Inheriting from this, in one and two dimension the localization length diverges in the limit w — > 
while in higher dimension the system is extended for sufficiently low w, see Sec. VA for further discussions.) 

Let us recall the canonical method of calculating the Green function - the path integral formalism [133] . Specifically, 



for wave propagation in an infinite random medium described by Eq. ( 14 ) , similar to quantum mechanics [1331 1134j , 



we may introduce the retarded (advanced) Green function G^i"^ [B] 



(u;^-i/)Gfl^(r,r') = 5(r-r'), 



(15) 



where ijJ± = lo ±i6 and 5 is a positive infinitesimal. The fluctuating dielectric e(r) follows Gaussian distribution with 



(6(r)e(r'))=A(r-r'), 



(16) 
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where A(r — r') = A(|r — r'|) is the correlation function. Then, one may cast the Green functions into the functional 
integral over some vector field, 0s(r), (f>*s{^)^ i-^-i 



/0..(r)C(rOe-^[^-ld[0:]d[0. 
J e~sl'Ps]d[(l,*]d[<f>,] 



where S[4>s] = =Fi / dnpKujj. — H)4>s is a bilinear action. (Throughout this review we use the notation or 
to denote the functional measure.) Here, the vector field (/)*(r), 0s(r) are either ordinary complex (with the 
subscript s = B) or anticommuting (with the subscript s = F) variables. The latter are called Grassmmanians. 
It should be stressed that for Grassmannians, the complex conjugate is purely formal: an anticommuting variable 
and its complex conjugate should be understood as independent variables. Notice that throughout this review the 
independent anticommuting degrees of freedom are even which allows us to freely move the measure under the integral: 
it can be placed either before or after the integrand. 

It is very important that the disorder, i.e., a;^e(r), enters into both the denominator (namely the normalization 
factor) and the numerator of the functional integral. This makes the subsequent disorder average a formidable task. 
Our aim therefore is to develop a first-principles theory such that disorder is eliminated from the normalization factor. 
This is accomplished by the supersymmetric trick |127|, I128j . Specifically, we promote 0s (r) to a two-component 
supervector (or graded) field, (/)(r), and 0*(r) to (/)^(r), the Hermitian conjugate of 0(r), i.e., (For simplicity, in this 
section we ignore the time-reversal symmetry of the Helmholtz equation (14).) 



(r) 



0F(r) 
Mr) 



0t(r)^(0*(r),0^(r)). 



(17) 



As the vector field 0F(r) {(psi^)) describes fermionic (bosonic) particles, (f>{r) represents a particle which is a mixture 
of fermion and boson: the term 'super' thereby follows. 

By using identity ( A31 1 in Appendix [A] we may rewrite the Green functions as 



±zGfJ^^(r,r') = 

Z[J\J] EE 



S^Z[J\J] 



SJF{r')S4ir) 



.7t,J=0 



(18) 



d[0^]d[<?!>] exp <^ - / dr (Ti(t>\r){Lul - i?)(/)(r) + J^(r)0(r) -I- 0^(r) J(r 



In the absence of the external source, = J = 0, the action is invariant under the unitary transformation in the 
supervector space: this is the so-called supersymmetry or Z2-grading. The most striking property of Eq. ( 18 1 is 
that Z[J^ = J = 0] is unity. Indeed, we have here a supermatrix, AI (cf. Eq. (A17l), as th e Gaussian kernel, with 
the matrix elements Mpp 



BB 



=F«(w± - H) and AfpB 



Mbf = 0. From Eq. (A30l, it follows that in the 



absence of sources, = J = 0, the partition function Z = sdetAf ~ det(aj^ — H)/dct{ujj. — H) = 1, where the 
denominator (numerator) results from the integral over the commuting (anticommuting) variables, 03,(/)b (0p,0f)- 
Here, 'sdet' stands for the superdeterminant (see Eq. ( A27) in Appendix|A]for the definition) and 'det' for the ordinary 
determinant. The striking property of Z[J^ — J = O] = 1 is due to equal number of anticommuting and commuting 
components. It is this property that (i) leads to a compact expression after disorder averaging (for disorders now 
enter only into the exponent of the numerator) , (ii) keeps the full effect of the disorder averaging, and (iii) renders 
the theory free of the analytic continuation problem which is encountered in the replica field theory [38l HH [89] . 

In fact, Z[J'^ = J = 0] = 1 is a special case of the general theorem discovered in Refs. [551 1135l[137j . According to 
this theorem, for a function f{x) satisfying /(-|-oo) < oo, we have 



d<j>u^f{.MA)) = /(o). 



(19) 



with ■\/(070) being the 'length' of the supervector (p (cf. Eq. (A16)). Remarkably, it states that (if the integrand is 



rotationally invariant in the supervector space,) the integral is given by the value of the integrand at the boundary of 
integration domain. To apply this theorem to present studies we expand (/)(r) in terms of the eigenmodes of H , i.e., 
4>{'r) = where (/fn(r) satisfies Hi^ni^) — a;^(/J„(r). Substituting the expansion into Z we find 



Z[,P = J = 0] = n / dc^lMnt 



= 1, 



(20) 



where the convergence of the integrand is guaranteed by the positive infinitesimal 5. 
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B. Origin of supermatrix field theory 



Above we explained the necessity of introducing functional integral over the supervector field if we would wish to 
get rid of disorders from the normalization factor. But this is not the end of the story. As will become clearer, it is the 
very beginning instead! Below we will adopt the heuristic approach of Bunder et. al. |138| to explain - in an intuitive 
manner - why eventually we will deal with an effective theory of supermatrix field, rather than the above-mentioned 
supervector field. To this end the source term is unimportant and we therefore set , J to zero. 

Let us perform the disorder averaging. As a result. 



Z 



{Z) = I d[0t]d[^]exp|- 



drdr' 



*0t(r')<5(r' - r) (V^ + u^) c^{v) + ^^A(r - r') (^^^^(r)) {ci^\v')c^{r')) 



The mathematical structure of the exponent is now undergoing a dramatic change: it is no longer quadratic in the 
supervector field. Rather, a quartic term appears. The latter may be viewed as the effective interaction among 
'elementary particles' represented by (j). 

On the other hand, experiences have shown that interesting physics in disordered systems arises from multiple 
scattering. This suggests that we should not treat the 'interaction' perturbtively. Instead, we have to keep track of 
its full effects. To this end, we recall the well-known identity: J dx5{x — xq) = 1. Suppose that its analog, 



i?Q(r,r')<5(Q(r,r') 



^) ^ (t>Hv')) = 1 



(22) 



(for all fixed r,r'), exists in supermathematics. Then, the Dirac function in this 'identity' enforces Q{y,y') to have 
the same structure as the dyadic product 0(r) (g) ^^(r'). Therefore, it must be a supermatrix, i.e., 



n(v rM= ( ^ff('''''') <9FB(r,r') 
' ^" WBF(r,r') gBB(r,r') 



(23) 



Here Qff,Qbb (Qfb,Qbf) are (anti)commuting variables. (More precisely, the former (latter) belongs to the subset 
{A^) of the Grassmann algebra and has even (odd) parity, cf. Appendix [a}) The measure is defined as DQ = 
Tr~^dQpBdQBFdQFFdQBB, where we have ignored the arguments r,r' to make the formula compact. 
Now, let us use the identity (A25) to rewrite Eq. (pTl) as 



{Z)= / d[0t]d[0]exp 



drdr'str 



i6{r' - r) (V^ 



- OJ 



q}\r') + — A(r - r')(0(r) ® <^^(r'))(<^(r') (b\v)) 



(24) 



where 'str' is the supertrace (see Eq. (A23) for the definition). Inserting the 'identity' (22) into it, we obtain 
(Z) = / d[0t]d[0] ( / D[Q]5{Q{r, r') ^ 0(r) ® cf,\r')) 



X exp 



drdr'str 



z(5(r' - r) (V^ + a;^) (j){r) (g) 0^(r') -I- ^A(r - r')(<?;'(r) (g 0^(r'))(0(r') ® 0^(r)) 



Suppose that it is legitimate to exchange the integral order. Integrating out the (/)-field first gives 
{Z) = I D[Q]J[Q]exp j / drrfr'str 



(25) 



z<5(r' - r) (V2 + u:^) Q{r, r') + - A(r - r')Q(r, r')g(r', r)J j , (26) 

where we have used the 'identity': / DQS{Q - (g (j}^)f{(j) ® (j)^) = J DQ5{Q - (p® <i)^)f{Q) and 

J[Q(r,r')] = j d4>\v')dcj,{v)5{Q{v,v')-^(v)®(l>\v')). (27) 

The latter is the analog of the well-known identity: / dx5{f{x)) = J^i \f'i^i)\^^ with f{xi) — 0. 

Thus, we have achieved an important result: disorder-averaged correlation functions can be traded to a functional 
integral over a supermatrix field. This is the core of the so-called superbosonization |138j . Nevertheless, the above 
derivations are intuitive and largely formal. In particular, we have paid no attention to the mathematical foundation of 
manipulations with the Dirac function of supermatrices. In fact, the superbosonization can be established on the level 



of mathematical rigor without involving the Dirac function of Grassmannians |138j . However, this requires advanced 
knowledge on supermathematics and we shall not proceed further here. In the next section, we will follow the more 
conventional route of using the super-Hubbard-Stratonovich (HS) transformation to introduce the supermatrix field, 
which was first done by Efetov |93) . 
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IV. SUPERSYMMETRIC FIELD THEORY OF LIGHT LOCALIZATION IN OPEN MEDIA 

In the remainder of this review, we wiU study wave propagation in large scales by using the supersymmetric field 
theoretic approach. Due to the interplay between wave interference and wave energy leakage through the air-medium 
interface, localization physics of classical waves in an open medium is even richer in comparison to that in an infinite 
medium. At the technical level, the corresponding field-theoretic description dramatically differs from that for infinite 
media. In this section, we will review the supersymmetric field theory of light localization in random open media 
developed in Ref. [ST]. We will first reproduce the supermatrix a model of Efetov in the context of classical waves 
in an infinite medium; Then, we will derive the air-medium coupling action that crucially constrains the supermatrix 
field on the interface. 

In reality, localization properties are probed by quantities such as the correlation function 3^(r,r';w), the wave 
intensity distribution [S3], the transmission distribution [55M38j . etc.. The microscopic expressions of these observables 
involve the product of the retarded and advanced Green functions. Thus, we have to double the supervector so as 
to account for the distinct analytic structures of these two Green functions: this defines the advanced/retarded ('ar') 



space, with index m = 1,2. Because Helmholtz equation (14) is invariant with respect to the time-reversal operation 
(the complex conjugation), we may consider the wave field and its complex conjugate as independent variables. To 
accommodate this degree of freedom we need to further double the supervector, which defines the time-reversal ('tr') 



space, with index t — 1,2. Finally, the doubling (17) defines the fermionic/bosonic ('fb') space, with index a = F, B. 



We have thereby introduced an 8-component supervector field ipir) = {V'mat(r)}. 

A. Nonlinear supermatrix a model for infinite media 

Let us start from the case of infinite media. It is necessary to reproduce the nonlinear supermatrix a model for 
light localization and discuss substantially the underlying technical ideas. The reasons are as follows, (i) Although 
the derivations for classical scalar waves are largely parallel to those of (spinless) de Broglie waves [51], below many 
detailed treatments are different from Ref. [88]. (ii) In doing so, we wish to show that the mathematical rigor of the 
nonlinear a model must be understood correctly. That is, it is not a 'rigorous' mapping of the microscopic equation 



(14); rather, it is an effective (low-energy) theory derived from the latter under some parametric conditions, (iii) 
Although the low-energy field theory turns out to be the same for classical and de Broglie waves, the conditions 
justifying such a theory are not. This leads to far-reaching consequences. In fact, as we will see below, localization 
physics of light differs from that of de Broglie waves in many important aspects. Therefore, it is necessary to have at 
our disposal a complete list of these conditions, (iv) Reflections of the derivations may provide significant insights into 
developing a first-principles localization theory incorporating vector wave character and (or) linear gain effects. Both 
are fundamental issues in studies of classical electromagnetic wave localization, (v) Although theoretical works using 
the supersymmetric field theory to study light localization have appeared, as mentioned in the introductory part, to 
the best of my knowledge, detailed derivations of this theory in the context of classical waves have so far been absent, 
(vi) We hope that a technical review on the nonlinear supermatrix a model of light localization in infinite media may 
help the readers unfamiliar with Efetov's theory to better appreciate how this technique unifies mathematical rigor 
and physical transparency, and powerful it is in studying light localization. 



Similar to Eq. (18), one may express the product of the retarded and advanced Green function in terms of certain 
derivative of the partition function, Z\j'^ , J], with respect to the source, J^(r), J(r). The partition function is now a 
hmctional integral over the supervector field, t/'j a-nd its Hermitian conjugate, V'^i 



Z{J\ J] = z J di^J^'Hi^] (exp |- J dv [iiP^K (V^ +uj^+ uj'^t{Y) - ww+A) + j'' K^j + i;^Kj] 



(28) 



where A = erf' (g) 1^ ig) 1*'', with af- , i = 0, 1,2,3 and =1'''^ the Pauli matrices defined on the space X ='ar', 
'fb', and 'tr', K is a metric tensor, and to make the formula compact we have dropped out the arguments of the 
fields. The numerical coefficient z depends on the metric tensor K. Similar to the example discussed in Sec. |III A[ 
this normalization factor is independent of disorder |139j . and its explicit value is not given here since it does not 
affect the subsequent analysis. We will postpone giving the explicit form of the metric tensor till Sec. |IV A 3[ At this 
mo ment we merely mention that the choice of K ensures the convergence of the above functional integral. In deriving 
Eq. (28) we have omitted all the cj^ terms, since we are interested in large-scale physics where the characteristic time 



(~ oiJj is much larger than the inverse (angular) frequency ui ^, i.e., 

w < w. (29) 
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To proceed further, we consider a simpler case for dielectric fluctuations, i.e., A(r — r') = A S{r — r') with A being 
the disorder strength. Performing the disorder averaging, we obtain 



exp <— dr 



ii;^ K {V^ + LU^ - uji:j+ 



A 4 -[ \ 



(30) 



Importantly, in the absence of the frequency and source terms, i.e., a)+ = = J = 0, the action namely the exponent 
is invariant under the gauge transformation U: tjj — > Utp, provided that U preserves the metric tensor, i.e., 

U^KU = K. 



Furthermore, by the construction of the supervector field it satisfies the 'reality condition', 

fb 



-ia*^ 



(31) 



(32) 



erf 

where the superscript, 'T', stands for the transpose. Its invariance under the gauge transformation requires 

= CU^C'^. (33) 



All the supcrmatrices U satisfying Eqs. (31) and (33) constitute the symmetry group of the system and is denoted as 
G. 



Effective interactions and low-momentum transfer channels 



Similar to Eq. (21), the action in Eq. (30) describes the dynamics of an 'elementary particle' which is a mixture 
of fermion and boson. Remarkably, the quartic term - arising from the disorder averaging - introduces the effective 
interaction between elementary particles. Written in terms of the components of the supervector, it is 

4 p 4 p 



2 



dr^, {-kK)l (^] ® 4'K^'J■) V,' , (34) 



with k = which accounts for the anticommutating relation of Grassmannians. Here, the subscripts ■ ■ ■ are 
the abbreviations of the index (mat) of the supervector components, and the Einstein summation convention applies 
to the index. Recall that if the elementary particle is coupled to an external field, (14xt)ij, then the action acquires 
an additional term J dr^plKij{Vey^t)jktpk (see, for example, the ujuj'^A term in the action). Comparing this structure 
with the interaction (|34|, we see that the particle self-generates some effective fields dictating the particle's motion. 
Specifically, the first equality suggests an effective scalar field ~ il^jKntfji (because this is a number), while the other 



two (which are related to each other via the reality condition (32|) of a supermatrix structure ~ ® K (because 
this is a dyadic product). 

In the Fourier representation, these self-generated fields are decomposed into slow and fast modes: the former 
fluctuates over a scale much larger than the mean free path /, with q = |q| (q being the Fourier wave number) much 
smaller than while the latter over a scale ^ I with q > . (For ballistic samples, where the mean free path 
and the system size are comparable, the fast-slow mode decomposition becomes very subtle, and it turns out that a 
low-energy theory is totally different from that to be derived below. We shall not discuss this issue further, and refer 
readers to the original papers ^143- 145^ and the monograph jl46j .) The slow mode possesses mathematical structure 
inheriting from the self-generated fields. Each slow-mode structure defines a specific low-momentum (i.e., ql <^ 1) 



transfer channel. To illustrate this we pass to the Fourier representation and rewrite Eq. (34 1 as 
^^■^ /■ j„ ^ ,t r.- ,^2 Aw* f dk dk' dq 

k+q 

k-l-q 



2 ./ (27r)<^ (27r)'' (27r)'* ^ V^^' k qy -/^ 

dk dk' dq , + / , , t tA , 
(2;0^(2^(2;0^^^>^^K'^'^W)^- 



{2'kY {2ttY (27r)" 
Aw* f dk dk' dq 
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k'+q 



-k+q I \-k'-q '''+1 * 

(a) (b) (c) 



(d) 



(e) 



FIG. 5: The (effective) particle-particle interaction may result in a small momentum (q) transfer via the direct (a), the diffuson 
or exchange (b), or the cooperon (c) channel. Successive small momentum transfer via the second (third) channel leads to 
a ladder diagram (d) ((e)), which is the prototype of the diffuson (cooperon) in the diagrammatic perturbation theory. The 
dashed line representing the interaction enforces two spatial coordinates involved to be identical, reflecting the 5-correlation 
nature of disorders. 



Here, tp\^ (■0k) is the Fourier transformation of i/j i"^^)} and stands for annihilating (creating) an elementary particle 
of momentum k (-k) . Notice that the interaction conserves the total momentum. Each equality above is respectively 
represented by an interaction vertex shown in Fig. [s] (a)-(c). The entity of two green (blue) lines defines a channel 
through which a momentum q is transferred. (With the momentum k or k' integrated out,) it varies over a scale 
^ l/q 3> I. Therefore, the entity specifies a slow-mode structure, and the interaction introduces the slow-mode 
coupling. 

Importantly, since the transferred momentum is small enough, i.e., ql 1, these three channels do not overlap. 
Following Ref. |133j we call the first, (a), the direct channel. Two particles undergoing scattering due to this interaction 
acquire a small momentum change q (or — q). The second, (b), may be called the diffuson (or exchange following the 
term of Ref. |133j ) channel and the third, (c), the cooperon channel. In fact, if the momentum q transfer through 
the second (third) channel occurs successively, then a diagram represented by (d) ((e)) results, which is the prototype 
of the well-known diffuson (cooperon) in the diagrammatic perturbation theory j741 1133] . Technically, the latter may 
be achieved by inserting Eq. (35) with corresponding low- momentum transfer channel into Eq. (30), expanding the 



quartic term and summing up the infinite series. 

The fast modes do not affect large-scale physics, and the interaction is dominated by the slow modes, 



dk dk' dq 

(2^(2^(2^ 



k+q, 



(36) 



with j' = Jgi^i^ where the factor two results from the reality condition (32) namely the reciprocal relation between 
the diffuson and cooperon channel. Define 



Jq^ 



dk 



v'0k^V'-k+q, 



dk 



k+q 



(37) 



we write Eq. ( 36 ) as 



(i^^K^j) drw 



4 



dq 



(27r 



(38) 



where in obtaining the second term we have used the identity ( A25 ) 



2. Introduction of the supermatrix Q-field 

The two exponents, Ii = exp ^f- /' JqJ-qj and I2 = exp |Au;* /' str(MqM_q)|, can be decoupled 
via the HS decoupling which transforms a quartic interaction into a quadratic one. The super-HS decoupling below 
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for the quartic interaction of the supervector field was introduced by Efetov (see Refs. 
and the original derivations). Here, we wish to adopt the treatments of Ref. |89j. 
By inserting the identity: 



l=J\flN, N= / I?[fq]exp 



1 



F F 



into l^ , we obtain 



/i =M- 



|i?[£Jexp|-^| ^£qf_,|exp 



2A 7 {I-kY 



for a technical review 



(39) 



J„J- 



(40) 



where f _q £*. Then, for the integral we make the change of variable: f q — ?► £q + iAw^ Jq which does not affect the 
convergence. As a result, 

= [ D[£^] exp r _ 1 /' -;5%7^qf-ql ■ (41) 



Therefore, we achieve the HS decoupling in the direct channel: on the right-hand side the action is linear in Jq. 
Conjugate to the 'external' field structure ~ ip''' Ktp which is a commuting variable (because Grassmannians appear 
in pairs), the decoupling field £q is a complex scalar field. This is the ordinary HS decoupling |133j . 
For the exponent I2, by inserting the unity: 



1^ I D[Q^]exp(^-A-' I ^str(gqQ_q)|, 



(42) 



where the supermatrix Qq has the same symmetry structure as the dyadic product -0 ip''^ K, and the convergence is 
assumed, we obtain 



l2 = J D[Q^] exp |-A-i J ^ str(gqQ_q) + Ao.^ J str(MqAf_q) 
Making the variable transformation: Qq — > Qq — Aw^Afq, we find 

l2= I D[Q^] exp |2c^2 J _dq_ st,(^.j^g_^) _ ^-1 J str(QqQ_, 



(43) 



(44) 



This is the super-HS decoupling in the diffuson-cooperon channel: on the right-hand side the action is linear in Mq. 
In contrast to Eq. (41 1, the super-HS decoupling in the diffuson-cooperon channel has a normalization factor of unity. 
This is a consequence of equal commuting (anticommuting) degrees of freedom which can be easily found to be 16 by 
explicitly writing down the matrix elements of the dyadic product -0 (8) i/;^ . 
Combining the HS decouplings ( 41 ) and ( 44 ) , we have 



AA-/i.[£q]eXp{-..^/'^/^£_q 



'-q'0k^'0-k+q 



2A 



—f f 



J i5[QJ exp 1-2^^1'^ J ^v4i^Q-q^-k+q-A-iy"^strQqQ_q 



(45) 



In the derivation above we have used the identity ( |A25 ). Most importantly, by the HS decoupling we transfer the 
interaction - a quartic term - into an action quadratic in tp. The price is the introduction of two new fields, the 
complex scalar field £q and the supermatrix field Qq. 



Let us substitute Eq. (|45| into Eq. (30). For the first line of Eq. (45), the first term in the exponent of the numerator 
is to locally renormalize the average refractive index, i.e., — >■ a?(l -t- £{r)), which is negligible. Therefore, the 
normalization factor f\f and the functional integral over £q cancel. That is, the decoupling in the direct channel gives 
a trivial factor of unity and plays no roles. As a result, upon passing to the real space representation, we obtain 



exp 



(46) 
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Notice that the Q-field is composed of slow modes and therefore varies over a scale much larger than the mean free 
path. Thanks to Eq. (31) the exponent of the right-hand side is invariant under the global (in the sense that the 
transformation is uniform in space) gauge transformation: ?/; — ^ U'ip, Q — > U^^QU, U E G. 

Then, we insert Eq. (46l into Eq. (30). For the action is quadratic in ip, we may integrate out this field exactly, 



obtaining (That to exchange the order of integration is legitimate is guaranteed by appropriate choice of the metric 
tensor K, see Sec. 



IV A 3 for further discussions. 



Z = I D[Q]t 



-fiq]-FjIq] 



Here, the action 



F[Q] = J drstr i^A-^Q^ - ^ In [-V^ - oj^ + ujCj+ A + 2iuj^ Q] 



(47) 



(48) 



Fj[Q] is the source-dependent action whose details rely on specific physical observable used to characterize the 
system's localization behavior and do not affect localization physics. Since we are interested in the general structure 
of the low-energy field theory in this section, we do not pay attention to its exp licit form. Importantly, although the 
functional measure formally includes both V'^ and ■0, the reality conditi on ( 32 ) reduces the indep end ent degrees of 
freedom by one half, which leads to the pre-logarithmic factor ^ in Eq. 
observables into an expression in terms of the functional integral over tl 



( 48 ) . Equations ( 47 ) and ( 48|) cast physical 
ic s upe rniatrix neld, Q(r) 



^ , . , , n the absence 

of the WW"*" A term and the external source namely Fj[Q], the action in Eq. (47) is invariant under the global gauge 
transformation: Q — >■ U~^QU where U G G is homogeneous in space. 



3. Nonlinear supermatrix a model 

We then consider the functional integral over the Q-fie\A. Again we ignore the source term, Fj[Q], for the same 
reasons as above. The program is first to find the saddle points (denoted as Qo) of the action F[Q]|cd+=o and then 
to integrate out (Gaussian) fluctuations around them. Let us substitute Q ~ Qq + SQ into the action -F[(3]|i:,+=o and 
expand it in terms of the deviation, 6Q. Demanding the linear term to vanish gives 

Qo{r)^^uj''Ago{r,r;Qo), go{r,v';Q) = i{r\{-W^ - lo'^ + 2tuj^Q)-y). (49) 

The imaginary part of the matrix Green function Qq implies that it decays exponentially when the distance |r — r'| 
is greater than a characteristic length which, as we will see below, is the elastic mean free path. In Appendix [B] we 
show that for low frequencies, i.e., 

ALJiy{u) < 1 ^ w < A'^/'^ (50) 



with i^(aj) ^ Lj'^ ^ being the density of states of free photons, diagonalizing Eq. ( 49 ) leads to two saddle points, qA and 
q{—kA), uniform in space (the 'mean field' saddle points), where q = jAujiy^ui). It is important to notice that these are 



not the complete solutions to Eq. (49). Indeed, this equation is invariant under the rotation, Qo ^ T ^QoT, T £ G, 



and the two saddle points, qA and q{—kA), are connected via this rotation. Therefore, the solutions to Eq. (49) 
constitute a manifold, 

Q^qT^^AT. (51) 

As shown in Ref. [89 , enforced by the convergence of the functional integral over ,Q and the analytic structure 
of the saddle point, the boson-boson (fermion-fermion) block in the fb-space is invariant under the non-compact 
(compact) group of transformations. This fixes the metric tensor (in its diagonal form). 



)fb 







K=[ Q ^"fb ) ® 1''-. (52) 



Combined with Eqs. ( |31[ ) and (33), this defines the symmetry group G = UOSP{2, 2\4:) for the orthogonal ensemble 
|139j . which is a pscudounitary supergroup. The notation to the left (right) of '|' stands for the metric in the bosonic 
(fermionic) sector: '2,2' refers to the hyperbolic metric, (-)-,-!-,—,—) and '4' to the Euclidean metric, (-I-, +,+,+). 
Importantly, there is a subgroup H G G which is a direct product of two unitary groups, UOSP(2\2), defined in the 
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advanced-advanced and retarded-retarded block, respectively. Elements of H generate rotation rendering A invariant. 
In this sense, elements of H may be considered to be identical. More precisely, in Eq. (51 1 T takes the value from the 
coset space, i.e., 



T e G/H = UOSP{2, 2|4)/C/05P(2|2) (g) UOSP{2\2). 



(53) 



As we will discuss in details in Sec. |IV C[ breaking of this global continuous symmetry leads to gapless collective 
modes, the so-called Goldstone modes. It is these modes, known as diffusons and cooperons in the perturbative 
diagrammatic technique [721 1731 176|. and their interactions that carry the full information on localization physics. 
This will be established in the following sections. We now proceed to analyze the action of these modes. 

Since we are interested in low-energy physics, the coupling between the Q-field spatial fluctuations and finite 
frequency (w+ ^ 0) effects is negligible, (The latter characterizes dynamic effects of fluctuations.) and they contribute 
separately to the effective action. To find the contribution of the former we set w+ = 0. Recall that the first order 
term in the (^Q-expansion vanishes. By keeping the expansion up to the second order we find the fluctuation action. 



{2ny 



str < uj' 



{2-kY 



^0 k 



5Q- 



(54) 



where we have passed to the Fourier representation: SQ{r) SQc 



GoiKQo)- (Notice that 
r', Qo)-) Suppose first that 



.q and Go (r, r', Qo) 

according to Eq. (49 1, Go (r, r', Qo) is translationally invariant, i.e.. Go (r, r', Qo) = ^^o (r - 
Qo = A. We see that fluctuations around it can be decomposed into two components: one commutes with it while 
the other anticommutes. Since the (anti) commutation relation is invariant under the global rotation, we may rotate 
A to arbitrary Qo in the saddle point manifolds. At the same time, we obtain two fluctuation components, the longi- 
tudinal (<5Qq) and transverse {SQV) components: the former (latter) commutes (anticommutes) with Qo. Physically, 



fluctuations along the transverse directions move Qo to somewhere the saddle point manifold (cf. Eqs. (51 1 and (53)), 
while fluctuations along the longitudinal directions bring Qo out of the saddle point manifold (see Appendix |C| for 
further explanations). 

Let us substitute the decomposition into Eq. (54). Integrating out the longitudinal components, (SQq, gives a factor 
of unity. We are left with a functional integral over the transverse component, 5Qq. Upon passing to real space, this 
gives the fluctuation action (see Appendix [C| for details). 



SFi = 



■nvDn 



rfrstr{V(T(r)-iAT(r))}2 



(55) 



Here, the transport mean free path determining the Boltzmann diffusion constant Do = l{u!)/d is of Rayleigh-type 
0162], i.e.. 



Therefore, we find that the inequality (50) is equivalent to the weak disorder limit, i.e., 

UjI{u!) > 1. 



(56) 



(57) 



It implies that the transport mean free path l{uj) is much larger than the wavelength X. Since we are interested in 
the large-scale physics, ui < Dq/P ~ the inequality (291 is guaranteed by (50). It is important to notice that 
here the weak disorder condition is established for low instead of high frequencies (w), in sharp contrast to electronic 
systems [55] . 



Next, we consider the contribution due to nonvanishing w. Keeping the w-expansion of the action (48 1 up to the 
first order, we find 



6F2 



drst 



rtr{A(r| (- 



J drstrjA^o (r, r; T-i(qA)r) } 



(58) 



where T is space-dependent. Ignoring spatial fluctuations of T and using the saddle point equation (49 1, we simplify 

it to 



(5F2 = i(wA)-iw+ j drstr{AT(r)-i(gA)T(r)}. 



(59) 
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Notice that F[Qo]\q+=o = 0. Finally, we reduce Eq. (47 1 to 
Z = J D[Q]e 



-fiq]-FjIq] 



Q(r) = r(r)-iAr(r), T(r) e G/H, 



(60) 



where the contributions to the action F[Q] (To avoid using too many symbols we use the same notation as Eq. (48).) 
are SFi 2, i.e., 



F[Q] = SFi + SF2 - J drstv{DoiVQf + 2iw+AQ}. 



(61) 



This is the nonlinear supermatrix a model action. 

Now let us briefly discuss the case where the time-reversal symmetry is broken (namely the unitary ensemble). In 
this case, the action (61 1 stays the same. The difference is the symmetry of the supermatrix Q. Indeed, the doubling 



introduced by time-reversal symmetry (in tr-space) is absent and the T-field is thereby a 4 x 4 supermatrix defined 
in the ar- and fb-spaces. Consequently, the reality condition (33) docs not play any role. (Correspondingly, the 



pre-logarithmic factor in the action (|48| disappears.) The system's symmetry is reduced, and Eq. (31) leads to the 
coset space G/H = C/(l, 1|2)/;7(1|1) (g) U{1\1). 



4- Ultraviolet divergence and long-ranged disorders 



Different from elastic waves, photons have arbitrarily small wavelength. Correspondingly, in high dimensions (d > 2) 
the right-hand side of the saddle point equation (49) suffers ultraviolet divergence (cf. Eq. (Bl ) and discussions on it 
in Appendix[B|. The divergence origins at that for large w the ^-correlated disorder is no longer a good approximation. 
Therefore, one may consider more realistic dielectric fluctuations namely long-ranged correlation of disorders described 
by Eq. (16]) |5|. Then, a fundamental issue of extreme importance arises: does this kind of disorders affect the 
universality of Anderson transition of light? This was first investigated by John and Stephen [5] by using the replica 
field theory. The key is to understand whether the low-energy field theory of localization is robust against finite-ranged 
correlation of disorders. Below, we will derive the low-energy supersymmetric field theory for this kind of disorders. 
(The corresponding replica field theory was derived in Ref. [5].) The program is very similar to that of deriving 
Eq. ( [61] ). Therefore, we will address the key differences and only present the final results. 

First of all, all the differences stem from the fact that upon performing the disorder averaging, the effective in- 
teraction is non-local which is cxp (to be compared with the 

quartic term in Eq. (30)). As such, the supermatrix field in the super-HS decoupling must be non-local also, since it 
has the same structure as the dyadic product, i.e., Q{r, r') ~ V(r) <8) ■tp\r')K. Consequently, the super-HS decoupling 
is modified to 



exp 



A(r - r') (V'^(r)ifV'(r)) {^Hr')Ki;{r')) dvdv' 



= j D[Q]expl^^ J (2cj2^t(i.)ifQ(r,r')^(r') + A-i(r-r')strQ(r,r')g(r,r'))drdr' 



(62) 



supermatrix field Q(v, r'), and its the saddle point equation is 

1 



Similar to the derivation of Eq. (48), with the ip field integrated out, we obtain a functional integral over the 

lu 

Qo(r, r') = ^oj^A{r - r')eo(r, r'; Qo), (63) 

where Qq has the same definition as that introduced in Eq. ( [49| . Here, (3o(r, i"') must be understood as the matrix 
elements of an operator non-diagonal in real space and correspondingly, Eq. ( 63 ) as an operator equation. To proceed 
further we introduce the Wigner transformation defined as Q(R, k) = f d{r — r')e^''' '^'"^'' ^Q(r,r'), with the center- 
of-mass coordinate R = (r -I- r')/2. Loosely speaking, R and k respectively describe the position and momentum of 
photons. Wigner transforming Eq. (63), we reduce it to 



A-i(*Vk)Qo(k) = ^c.2g„(k; Qo(k)), 



(64) 



where GoO^', Qo(k)) is the Wigner transformation of Qo. This leads to a diagonal mean-field saddle point, q{k)A, and 
modifies the saddle point manifold to 



Q(k) = q{k)T-^AT, T e G/H. 



(65) 
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Here, T depends neither on k nor on the center-of-mass coordinate R, and q{k) is a real positive function of k. 

As before the fluctuation action results from the transverse fluctuations, i5(5*(R, k) = [Qo(k), 5T(R, k)], around the 
saddle point Qo(k), where (5T(R, k) parametrizes fluctuations in the coset space G/H. The fluctuation action is 



dk 



{2ttY J (27r 



str{(5Q^(k)7i^JQ*_q(k)} 



(66) 



with 5(5q(k) being the Fourier transformation of 5Q*(R, k) with respect to R. Here is a 'Hamiltonian' depending 
on the 'external parameter' q. In the k-representation, it reads 



?i^(k) = A-i(zVk) + Lu^Go (k + |; Qo(k)) (k - |; -Qo(k) 



(67) 



where the first (second) term may be viewed as generalized kinetic energy (potential). 

Different from the case of (5-correlated disorders (cf. Appendix [C|) , the integrand in Eq. (66) does not vanish at 
q — > in general. Instead, (5(5* (R, k) are composed of low- and high- lying modes. To single out the former we have 
to expand 6Ql,(k) in terms of the eigenfunctions of the Hamiltonian Hq(k) denoted as w/(k; q) which satisfies 



H*(k)«i(k;q) = A,(q)«Kk;q) 



(68) 
for the 



with A/(q) being the eigenvalues. Here the subscript I — stands for the 'ground state' and I = 1,2 
'excited states'. Following John and Stephen 5j one finds that for general expression of A(x), at q = Eq. (64) 



leads to a unique ground state tio(k;0) oc q{k) with zero angular momentum. The ground state eigenvalue vanishes 
as Ao(q) = C2q^ hi the limit: q — > 0, and all the excitation modes A;>i(q) are gapped. Here, the coefficient C2 may 
be found by applying the Raileigh-Schrodinger perturbation theory to Eq. (68). Therefore, the low -lying transverse 
components must have the general form, (5(5*(R, k) — [Qo(k), <5T(R)]. Substituting it into Eq. (66) and performing 
the hydrodynamic expansion, we obtain 



dFi 



Trh'{uj)Do{uj) 



where the density of states is given by 



2cj, 



-Im 



dRstr{V(T(R)-iAr(R))}2 



dk 1 

{2Try k^ -LU^+ 2iuj^q{k) 



(69) 



(70) 



and the (bare) conductivity i^(a;)£)o(w) cx C2. Equations (|64| and (70) lead to the self-consistent Born approximation 



(or the 'coherent potential approximation' called by John and co-workers [11[5]) of the density of states 

The procedure of deriving the frequency contribution, dF2, is similar to that of deriving the second term of Eq. ( 61 ) 
Now, Eq. (58) becomes 



6F2 



tUJLJ 



str{Ago (k;T-i(R)(g(k)A)r(R))} 



dR 



dk 

{2'kY 
dk 

^str {AT-i(R)g„ (k; g(k)A) r(R)} 
dRstr{AT-i(R)Ar(R)} , 



(71) 



where in deriving the last line we have used Eq. (70). From Eqs. (69) and (71) we conclude that the nonlinear 



supermatrix a model action (61) is valid even for finite-ranged disorders, provided that the Hamiltonian "H^ has a 
unique isolated ground state. Tlotice that the explicit expression of Dq{uj) is generally changed. 



5. Possible generalizations 



Several exact solutions for strong localization in weakly disordered wires have been obtained by using the nonlinear 
a model (61), notably, the correlation function defined in Eq. ([!]) in QID infinite wires [MIj the first two moment of 
conductance ^7J HH] and transmission distribution [351 [36, 98 in finite wires coupled to ideal leads. These results 
have revealed a deep connection between the noncompactness inherited from the hyperbolic metric in the bosonic 
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sector and strong localization. A further example is the anomalously localized states in diffusive open samples [140]. 
As mentioned above, the hyperbolic nature of the metric stems from (i) the convergence of functional integral and 
(ii) the analytic structure of saddle points. To build a localization theory incorporating linear gain effects [53], one 
must inevitably treat both (i) and (ii), since now the sign of the convergence generating factor, S, is reversed, i.e., 
w+ — 7> uj~ . As such, discussions parallel to those of Ref. [53] must be made cautiously, and modified theory may help 
to study many interesting phenomena such as the gain-absorption duality |14H I142j . 

It seems possible to generalize the present formalism to work out a first-principles theory for electromagnetic wave 
(vector wave) localization, a long-standing issue in studies of light localization. Reflecting the scheme outlined in this 
section, one may naturally enlarge the supervector field V to accommodate the vector field structure. However, this 
step - right at the beginning of the entire field theory formulation - is by no means trivial. In fact, the additional 
field components are not trivially independent, and a crucial step would be to embed the constraint reflecting the 
transverse field character of electromagnetic waves into the functional integral formalism. Then, one may derive a 
low-energy field theory by following the above scheme. 



B. Field theory for open media 



We have so far focused on infinite (bulk) media, where the low-energy theory (61 1 is translationally invariant. The 
presence of the air- medium interfaces breaks this symmetry. In this part, we will show that the Q-field is constrained 
on the interface. The low-energy field theory thereby obtained describes many exotic localization phenomena, which 



will be exemplified in Sec. [v| and VI In fact, open meseoscopic electron systems (e.g., quantum disordered wires or 



small quantum dots coupled to ideal leads, superconducting-normetal hybridized structures, etc.) have been studied 
extensively (see Refs. [551 IMl [T^ for reviews). The effective field theory describing these systems has been worked 
out in Refs. [551 [551 [1471 1148j . All these works study physical observables such as the conductance or transmission. 
Whether localized open systems may exhibit macroscopic diffusion (a kind of hydrodynamic behavior) was not studied 
in these works. In addition, the systems considered there are in quasi one (e.g., quantum wires) or zero (e.g., quantum 
dots) dimensions. 

In classical wave systems, the experimental setup is generally more complicated. For example, it may be in 
higher dimensions [T^ 1149) . In addition, internal reflection may exist on the air-medium interface due to the 
refractive index mismatch [32]. The subjects to be investigated are more diverse. Besides of the physical observables 
mentioned above, particular attention has been paid to many other issues, e.g., time-resolved transmissions, coherent 
backscattering lineshape, and spatial resolution of localized modes. In Ref. [17] j the method developed for mesoscopic 
electronic devices [88, 98, 147, 148 was first generalized to open classical wave systems to address a number of issues 
of classical wave propagation in random open media. 



1. The interface action 

The propagation of classical waves in the random medium (the space supporting which is denoted as V_|_) is described 
by the effective Green functions, G^i (r,r'), r,r' € V+. These two Green functions differ from G^i (r,r') only in 



the domain of r,r'. According to the definition (15), for G^i'^{r ,r') the variables r, r' can be either in the air (The 
corresponding space is denoted as V_.) or the random medium V+. We also introduce an auxiliary Green functions, 
g^'i (r , r'), defined in the air satisfying 

(4-^)5y(r,r')=<5(r-r'), fff^-^(r, r')|rorr'ec = 0, r,r'eV_, (72) 

where C is the air-medium interface. To proceed further, we need to introduce a theorem due to Zirnbauer [148] 
and refined by Efetov [55], originally established for describing the coupling between leads and mesoscopic electronic 
devices. In the present context of classical wave systems, the theorem may be stated as follows: 

(Zirnbauer-Efetov) The effective Green function Q^i'^{v ,v') (r,r'E V^) solves 

(..2 „H±zBJc) ^y(r,r')-<5(r-r'), V„(r)ey(r,r')L,^ - 0, reC. (73) 
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Here, n(r) is the normal unit vector on the interface, pointing to V+, and the operator 13 is defined as 

{Bf){v)^ ( dr'Im[i?(r,r')]/(r'), i?(r, r') = V„(r) • V„(r-)g5 (r, r') , (74) 

JC 

with / a test function on W'-. Finally, the operator 5c is defined as / drSc /(r) = ^^^^c ^''-^/('"-l) ■ 

The th eore m states that wave propagation in random media follows a modified microscopic equation. It differs 



from Eq. (15) in that the effective 'Hamiltonian' acquires a correction, the iBSc term which is purely imaginary and 
located on the interface. The correction accounts for the wave energy leakage through the interface and therefore for 
the dissipative nature of quasi-normal modes. Correspondingly, by repeating the procedures of Sec. IV A 2[ (We refer 



readers to Appendix pP] for details.) we find that the action (61) acquires a correction, the interface action Fint[Q] 
given by 

where the internal reflection coefficient i?k^(r) generally depends on the transverse momentum, kj^. For simplicity, 
we ignore the k_L-dependence of the internal reflection coefficient from now on, i.e., i?kj_ = Ro- In other words, we 
shall not study the effect arising from that internal reflection depends on the angle which the k- vector makes with the 
interface. Because both i?o(r) a-nd Q{r) vary over a scale larger than I, in the limit ujI ^ 1 one may further expand 
the logarithm and simplify the interface action to [67j 



where (1 — i?o(r))/(l + i?o(r)) characterizes the strength of the air-medium 'coupling'. (This interface action (76) 
differs slightly from that derived in Rcf. [57| because internal reflection is small there.) 



2. Boundary condition satisfied by the Q-field 

To proceed further we use the so-called boundary Ward identity [126J. It states that the functional integral 
representation of the local observable, / £'[Q](- • • )e^^['31^^'"*['31, with (• • • ) depending on specific observables whose 
details are unimportant for further discussions, must be invariant under the infinitesimal boundary rotation: Q — >■ 
e~^^Qe^^ fa Q — [6R,Q]. Therefore, it is required that the resulting variation of the interface action identically 
vanishes, 

5Fi„t = y'dr,5cstr|<5i?(^^QV„(r)g-H^^^[Q,A])|^0. (77) 
For SR is arbitrary, it gives [150] 

(2cgv„(r)g + [g,A]) |^^^ = o, (78) 



where C is given by 



2 1 - i?n N, 



d-l 



Equation ( 79 ) reproduces the well-known result for the extrapolation length obtained by completely different methods 
[121H 123'. However, it should be stressed that the previous methods are valid only for diffusive samples, while the 
present theory is valid for both diffusive and localized samples. In the absence of internal reflection, i.e., Rq = 0, 



Eq. (79) gives the value of | ^ in three dimensions, in agreement with the result of Refs. [Ml I104j . This value is 
closeato the exact value 0.71/, obtained by solving the radiative transfer equation |103[ llSlI I152j . (This small 
deviation is due to the simplification of R\^^ = Rq made above. The exact value can be reproduced by the present 
field theory by retrieving the kj^-dependence of the internal reflection coefficient.) In a study of one-dimensional 
iS'S"5'-superconductor structures |153j . Kuprianov and Lukichev derived a boundary condition for the Usadel equation 



[1541 1155j similar to Eq. ( 78 ) , in which the matrix Green function acts like the supermatrix Q-field. 



Equation (78) is the boundary condition satisfied by the Q-field. Together with the nonlinear a model action (61) 



it lays down a foundation for quantitative analysis of wave localization in open media. The boundary constraint 



24 



plays an essential role in a field-theoretic description of localization in the presence of the interplay between wave 
interference and wave energy leakage. It is valid no matter whether the sample is diffusive or localized. Noticing that 



QV„(r)(5 mimics the boundary current, we find that Eq. (78) resembles the so-called 'radiative boundary condition'. 
The latter has been used by many authors to implement the normal diffusion equation to study light propagation in 
diffusive samples [12111122] , and was justified in the context of one-dimensional random walks |156j . Indeed, one can 



show that it is a special case of the more general constraint ( 78 1 



C. Spontaneous symmetry breaking and Goldstone modes 



The supermatrix model ( 48 ) gives rise to spontaneous symmetry breaking, and the reduced action ( 61 ) characterizes 
the energy arising from proliferation of the low-lying Goldstone modes associated with broken continuous symmetries. 
To see this let us first recall a more familiar model namely the Heisenberg model of ferromagnetism [15711158] described 
by the following Hamiltonian (in unit of temperature) , 

Jf[S] = -J^Si-Sj-h-^S„ S^ = l. (80) 

Here without loss of generality the magnitude of the spin S is constrained to unity, the coupling constant J > 0, 
and h is the external magnetic field. In the absence of h the system is invariant under the global three-dimensional 
rotation giving the symmetry group G — 0(3). Below the critical temperature, this symmetry is broken and all the 
spins point to the same direction, say the north pole [158j . Hence the adjective 'spontaneous' follows. Although 
in general individual ground state (at h = 0) breaks the symmetry, all the ground states constitute a degenerate 
manifold and restores the full symmetry G. Now, an infinitesimal magnetic field h — >■ plays a remarkable role: it 
selects a particular ground state from this manifold. In the selected state all the spins are parallel to h, since globally 
rotating this state costs infinite energy in the thermodynamic limit. This is the so-called spontaneous symmetry 
breaking, see Fig. |6] (a) . Suppose that the ground state undergoes local rotation and the resulting spin configuration 
varies over a spatial scale ^ . The transformed state must have a larger energy. The energy increase diminishes 
as the rotation becomes uniform, i.e., q 0. That is, spatial fluctuations of spins generate soft modes in the sense 
that their energy vanishes as g'' at g — )■ where the exponent v > 0, i.e., these soft modes are gapless. These 
gapless modes (sketched in Fig. [6] (b)) are the so-called Goldstone modes. However, it is important to notice that 
not all the inhomogeneous rotation costs energy. In fact, the ground state is invariant if the local transformation is 
generated by the two-dimensional group of rotation, H = 0(2), a subgroup of the full symmetry group G — 0(3). 
Therefore, the Goldstone modes are associated with the breaking of the continuous symmetry belonging to the coset 
space G/H — 0(3)/0(2) instead of the full symmetry group G (cf. Fig. |6](c)). 

The supermatrix a model bears a formal analogy to the Heisenberg ferromagnet. For the action ( |48[ ), in the absence 
of the frequency term ~ woj+strAQ it is invariant under the global gauge transformation: Q — U^^QU where [/ is a 
constant matrix belonging to the symmetry group G, which is U0SP{2, 2|4) for GOE systems and [/(1, 1|2) for GUE 



systems. Each 'ground state', Qg, solving the saddle point equation (49) breaks this continuous symmetry: in general. 



U ^QoU Qq. Now, if we restore the frequency term, i.e., cjcD+A with cD — >■ in the matrix Green function Qo (cf 



Eq. (49)), the system has a unique 'ground state' oc A. In other words, the frequency (or the convergence-generating 
factor S a uj — 0) term plays the role of the 'external field' and causes the spontaneous symmetry breaking. 

Similar to Heisenberg ferromagnet this 'ground state' is invariant under the local gauge transformation: Qq — > 
U-\r)QoU{r) where U{r) e H with H = UOSP{2\2) ® C/OS'P(2|2) for GOE systems and H = C7(l|l) ® ;7(1|1) for 



GUE systems. Then, T(r) G G/H describes the Goldstone modes and the action (61) describes the energy required 
for exciting low-lying modes. As expected -F[Q]|(i+=o vanishes if Q is homogeneous. 

V. MACROSCOPICS OF WAVE PROPAGATION IN OPEN MEDIA: GENERAL THEORY 

In Ref. [159[ . it was for the first time found at the full microscopic level that the diffusion coefficient in open 
media is position-dependent. However, very specific random media were studied in that work. Later, the microscopic 



approach used there was generalized to realistic random open media (see Sec. IV I and was used to systematically 
study transport of classical waves in these systems ^TUB^j. Armed with the first-principles theory developed, it was 
found that the wave intensity correlation - a key characteristic of wave propagation - in open media is very different 
from that in infinite media [72j [73l [94j. Specifically, such correlation function was found to obey Eq. ([4| in the 
frequency (lu) domain. Most strikingly, in Ref. [69[ it was shown analytically and confirmed numerically that in the 
one-dimensional case, the static local (position-dependent) diffusion coefficient exhibits a novel scaling behavior. In 
this section, we shall establish the scaling theory in the more general case of slab open media in high dimensions. 
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FIG. 6: Spontaneous symmetry breaking in a Heisenberg ferromagnet. The ground state breaks the global rotational symmetry 
G — 0(3) (a). Collective fluctuations around this state, the Goldstone modes, are gapless (b). They are associated with the 
breaking of the symmetry G/H = 0{3)/0{2) ^ (c). 



From the technical viewpoint, the source of the difficulty can be seen in that even near the critical dimension the 
ordinary renormalization group scheme does not apply: one cannot replace the Boltzmann diffusion coefficient by 
the renormalized one obtained for infinite media because it does not account for effects arising from the openness 
of media such as resonances [iBHSU] , We therefore resort to calculating the wave intensity correlation by using 
the developed first-principles theory and then proceed to establish the macroscopic (or 'hydrodynamic', the term 
commonly adopted by condensed matter physicists) diffusion equation describing wave propagation. Restoring the 
source term and taking the corresponding derivatives of the partition function (47), we cast the correlation function 
^ into a functional integral over the supermatrix field Q(r), i.e., 



8w 



D[Q]stT{A+Q{Y)A^Q{Y'))e 



-F[Q] 



(81) 



for (I) ^ tj, where A± 
open media. 



(T3 (g) (1 ± cr|'') (g) (1 — 0-3''). This expression is valid for both infinite media and finite-sized 



A. Wave propagation in infinite media 



It is very illustrative to calculate Eq. (81) first for infinite media. In doing so, one may become familiar with 
perturbative calculations of the supersymmetric field theory, and appreciate - at the perturbative level - the main 
technical differences between infinite and finite-sized open media. We shall explain in detail the calculations for GOE 
systems and give the results for GUE systems directly. 



Notice the following identity. 



1. Parametnzation of the supermatrix t 



T . 



Til T12 
T21 T22 



Til 
T22 



1 Tfi^Tia 
T22 T21 1 



Til 
T22 



(l + iW). 



(82) 
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FIG. 7: In the perturbation theory of the nonhnear supermatrix a model, each diagram is composed of the bare propagator 
(top) and the interaction vertexes (middle and bottom). 



Upon inserting it into Q, we find 



{l + iW)-'A{l + iW), W=(^l^^] , (83) 



which is the so-called rational parametrization. The structure imp osed on the W matrix is fixed by two conditions. 
That is, W anticommutes with A and = KQK due to Eq. (31). The 4x4 matrix B can be parameterized as 

where the entries of a's and 6's (a's and 77's) are complex commuting (anticommuting) variables. They are not 
independent due to the constraint. 



introduced by the condition (|33| , and the number of independent integral variables is reduced by one half. In Eq. (|84j) , 
the subscript i = 0,3 (1,2) stands for the diffuson (cooperon). An advantage of the rational parametrization is that 
the corresponding Jacobian (the so-called 'Berezinian') is unity. For perturbative calculations below, it is sufficient 
to take this as an exact mathematical theorem. We refer the reader interested in its proof to Ref. [S8 . 



2. Wave interference corrections 



For large w we may expand Eq. (81 1 in terms of W. The action F[Q] is given by X^^^i [W^] with F„[Vl^] oc VK^". 
The pre-exponential factor leads to another VF-expansion, J2Ts=a^^''^i^+-^r[W{r)]A^Pg[W{r')]), where — 
0(\V). As a result, 

3^(r,r';^) = ( V str(A+P.[W^(r)]A_n[W^(r')]) V £y I VF„[T4^] I ) , (86) 



8a; 



/ 00 00 /00 \ \ 

J2 striA+Pr[W{r)]A^Pg[W{r')]) ^ - ^F„[T4^] ) . 

\r,s=0 m=0 ' \n=2 / / q 



where the average is defined as {■}o = J D[W]{-) exp{-Fi[W]), with Fi[W] = - f J dratT{Do{\/iW)^-iu){iWf). This 
perturbative expansion describes the interaction between diffusons and cooperons, and each term can be calculated 
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(a) 



(b) 



FIG. 8: Diagrams give the leading order wave interference correction to the Boltzmann diffusion constant. 



(87) 



by using Wick's theorem with the help of the following contraction rules, 

AT:v{W{r)MW{v'))Q = yo{r,v';Co) [(strM + M) - A(str(MA) + MK)] , 
47ri/(str(Ty(r)M)str(H/(r')M'))o = 3^o(r, r'; w)str((M - M){M' - AM' A)), 

with M — KCM"^ C'^ K . The propagator, 3^o(r, r'; tl'), solves the normal diffusion equation, 

{~iCo - V2)3^o(r, r'; Hj) = 5{v - r'), 

which, importantly, is translationally invariant, i.e., 3^o(r, r';w) = 3^o(r — r',0;a)). The contraction rules can be 
proved straightforwardly by using the parametrization (83) and (|84|. With all the contractions made, each term 



in the expansion (86) can be represented by a diagram composed of the (bare) propagator 3^o(r,r';a;) (representing 
the corresponding contraction), and the 'interaction' vertices -Fn>2[W^] and Pr[W^] (Fig. 7]). The diagram is either 
reducible or irreducible, depending on whether it may be divided into two disconnected subdiagrams by cutting some 
bare propagator. 



Let us start from the lowest order term {0{W )) in Eq. (86), i.e., 

2 



■Kv{uj) 
4a; 



D[W]e~^^^'^hiT:{A+iW{Y)A^iW{v')). 



(89) 



With the help of the contraction rules ( |87[ ) we find that at this perturbation level, the correlator 3^(r, r'; Cj) is identical 
to the bare correlator 3^o(r,r';w), which excludes all wave interference effects. The subleading term is OiVV^) with 
all the irreducible diagrams given in Fig. [sj (Notice that Fn[W] oc VT^" = ©(^^^(n-i)-) j^j. ^^j^j ^ Q^^y^^ rpj^^g^ 
terms contribute wave interference corrections to the bare correlator. More precisely, at this perturbation level the 
correlator is (see Appendix |e] for details) 



yQ{r,r'-Qj) + 5yi{v,v'-[b) 



iii'{u)) 



where 



(90) 



(91) 



By re-exponentiating the last factor (which is equivalent to summing up all the reducible diagrams with the one-loop 
vertex) we rewrite Eq. ( 90 ) as 



3^(r,r';c:;) 



D[W] 



J dristi(D(uj)(ViWf-ii:i(iWf) 



s\:T:{A+iW{v)A^iW{Y')). 



(92) 



Comparing the Gaussian wight with that of Eq. ( 89 1 , we find that wave interference renormalizes the bare diffusion 
coefficient, i.e., Dq -5> D{Cj) = Dq + 5D^^\Cj), and 5D'^^^{Cj) is the well-known weak localization correction [72l [73l [76] 
namely Eq. (|3|. The correlator ( [92| ) solves the macroscopic equation ([2|. 

Importantly, the diagram (b) in Fig. [s] guarantees that the wave interference correction to the bare propagator 
3^0 vanishes when J^o is homogeneous in space [H7] (see also Appendix [e] for the proof). In general, the roles of the 
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expansion of the pre-exponential factor are to conform with an exact Ward identity reflecting the energy conservation 
law, 



— lUJ 



(93) 



This must be respected at all levels of the perturbation. In the diagrammatic perturbation theory, conforming to this 
identity necessitates the introduction of the well-known Hikami box |102j . 

For GUE systems the leading order weak localization correction (91) vanishes. This certainly does not imply the 
absence of interference effects in the weak (not strong) localization regime. The point is that wave interference (see 
Fig. 11 for example) can be much more complicated and quantitatively they give higher order weak localization 
corrections. It must be emphasized that these higher order corrections cannot be obtained by calculating reducible 
diagrams with one-loop vertex. 

Let us now consider the terms ~ 0{W^). The irreducible diagrams of this order are given in Fig. [9] (Recall that 
F2[W] = 0{W^) and F^lW] = 0{W'^).) The derivation is similar to that of deriving 6D'^^'>{Cj) but is much more 
complicated. By using the contraction rules (|87| one may show that for GOE systems the correction given by Fig. [9] 
(denoted as 5y2{r,r';uj)) vanish, i.e.. 



SD^'^\Cj) = 0, for GOE. 



(94) 



For GUE systems here we will emphasize the key points and only give the final result |160j . First of all. Fig. [9] (d) 
(g) again conforms to the Ward identity (93|. Secondly, dy2{r,r';uj) consists of two contributions, i.e.. 



^y2(r,r';w) 



TTh'{iu)\ nvDo 

1 ( 3^o(ri,r2;w) 
2 



D[W]e-^^^^^siv{A+iW{r)A_iW{Y')) j drij drzstr {^"iW{ri)V^HW{v2)) 



5{vi - V2)5c.p + i?oV"V^:Vo(ri, r2; c^) 



/ yo(ri,r2;tD) 

\ TTV 



(95) 



The first contribution leads to a wave interference correction to Dq which is local in space thanks to the factor 
5(ri — Y2), while the second one to a correction non-local in space. The non-locality arises from diagrams such as 
Fig. [9] (a). Because the system is isotropic, the crossing terms, ~ V^V^, a 7^ /3, do not contribute upon integrating 
out the space coordinates. As a result, Eq. (95) is simplified to 



5y2{r,v'-i::) = - 



1 /^3^o(ri,r2;(I'' 
2 



D[W]e-^^^^^stT{A+iW{r)A_iW{r')) J dv^j dr 2Siv\/ iW [r i) -^/iWiv 2) 



5{vi - V2) 



Do 



V^3^o(i-i,r2;w 



/3^o(ri,r2;(I') 



(96) 



In the hydrodynamic limit the spatial and time derivatives in the macroscopic equation are decoupled. Taking this 
into account, Eq. (96) reduces to 



53^2 (r,r'; (I;) 



1 / yo(ri,r2;a)) 

2 V TTV 

2 



D[W]t 



-FilW] 



str{A+iW{r)A^iW{r')) J dvij dr2StrViW^(ri) • \/iW{r2) 



6{ri - ra) -I- -{DoV^ + iCb)yQ{ri,r2;^) 



3^o(ri,i'2;w) 



) / D[W]e-^'^^^sti{A+iW{r)A^tW{r'))'^ f dri6D^^\(Jj)sir{\/tW{ri)f 



(97) 



where the wave interference correction to Dq is local, read 



1 1 

2 ~ d 



^0 



po(o,0;^) 



for GUE. 



(98) 



At this perturbation level, the correlation function is 3^ = 3^o + <53^i + (53^2 which solves the macroscopic diffusion 
equation ([2]), with the dynamic diffusion coefficient replaced by £'(aj) = 1^0-1- 8D'^'^\lS) + 6D'-^\u}). At the higher 
order perturbation level, Eq. ([2]) remains valid with I?(a)) including corresponding higher order corrections. 
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(a) 



(b) 



(c) 








(d) 



(e) 



(f) 



(g) 



FIG. 9: Diagrams give the subleading order wave interference correction to the Boltzmann diffusion constant. 



Notice that there are two opposite contributions to (5D*^^^(w), accounting for the factor ^ and — respectively. The 
wave interference picture underlying the former |161l 1162] (see also Fig. 11 below) may be related to the interesting 



phenomenon of coherent forward scattering }163j . Importantly, the factor — ^ ) ensures the renormalizability: ap- 
proaching two dimensions (from below), this infinitesimal factor cancels one logarithmically diverging factor arising 
from 3^o(0,0;u;) and makes the two-loop correction diverge only logarithmically. Upon making the analytic continua- 
tion from dimensions below two to above, we find from Eqs. (91 1, (94 1 and (98) the well-known scaling equation for the 
dimensionless Thouless conductance g{Lcj) — tti'DqL'^''^ with Lcj = Dq/lj and D,d the volume of d-dimensional 
unit sphere, i.e.. 



ding 
dlnLr, 



for GOE, 



d~2 



1 

2g2 



O(^), forGUE 



(99) 



near the critical dimension [451 

We wish to make some remarks here, (i) In the presence of magneto-optical effects the one-loop weak localization 
may be fully suppressed. In this case, the critical theory is described by the second equation (99). This overcomes 



the notable difficulty of the self-consistent diagrammatic approach |75J [73] which cannot be used to study localization 
of GUE systems. This is a simple example showing the field theoretic approach to be a powerful technique in studies 
of localization in magneto-optical systems, (ii) The criticality described by the scaling equation is quantitatively 
accurate in2-fe (e = d — 2 which should not be confused with the dielectric constant) dimensions since the critical 
point g* = 0(6"^) > 1. One may extrapolate the results into three dimensions by setting e = 1. In particular, 
we find that the loffe-Regel criterion is also valid for GUE systems. This may guide experimental search of light 
localization in magneto-optical systems, (iii) In combination with Eq. (56) this criterion predicts that there is a 
mobility edge separating low-frequency extended states and high-frequency localized states, (iv) In low dimensions 
(d < 2), the dimensionless conductance decreases as the frequency uj is lowered. In this case, all the higher order 
irreducible diagrams are important. Similar to the above analysis, these give local corrections to the Boltzmann 
constant. This implies the validity of the macroscopic equation ^ for low Co. However, obtaining the explicit form of 
D{uj) is beyond the perturbation scheme and one must resort to non-perturbative treatments of the supersymmetric 
nonlinear a model. This indeed has been done by Efetov and Larkin in the QID case, (v) According to the 
exact (non-perturbative) solution given in Ref. 194j, in QID the localization length f ~ v{uj)DQ{io). This implies 
that the localization length diverges as in the limit w — > 0. Notice that this scaling law is universal (but the 
coefficient not): it is the same for both GOE and GUE systems, (vi) In two dimension the localization length diverges 



exponentially in the low-frequency limit. The behavior depends on system's symmetry: In^ 



and In ^ - 



for GUE systems. The latter results from the infinitesimal factor of (^ — ^). 



for GOE systems 



B. Wave propagation in open media 



We now consider the case of open media and investigate how the results above arc modified. Let us begin with 
a pictorial illustration in terms of the supermatrix field of the main difference between infinite and finite-sized open 
media. We then substantiate the picture with explicit calculations of the wave intensity correlator. 
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TABLE I: Ferromagnet-supermatrix cr model analogy 





Heisenberg ferromagnet 


supcrmatrix a model 






F\0] 


spin 






constraint 




Qi = 1 


coupling constant 


J 




magnetic field 


h 


4 


ordered phase 


ferromagnet 


metal 


disordered phase 


paramagnet 


insulator 



1. A ptctonal illustration 

For simplicity let us consider the discrete space lattice (with the lattice constant set to unity). In this case, the 



nonlinear a model action (61 ) is written as 



F[Q] = -i^Do J2 ^^^QrQ] + ^^^^ strAQ,, Qi = 1, (100) 

(ij) » 



where i stands for the lattice site and the summation is over all the nearest neighbors, {ij). Comparing Eqs. (80) 
and ( 100 ) we find that these two models bear a firm analogy summarized in Table [ij Notice that the supermatrix 



structure of Q as well as the supertrace operation is irrelevant for present discussions. 

Identifying this analogy, we may translate some well-known facts regarding the ferromagnetic system |158) to the 
present supcrmatrix a model. For the ferromagnetic system, if the temperature is sufficiently low, spins tend to align 
in parallel as a result of the so-called spontaneous symmetry breaking. In the presence of the external magnetic field 
h they are all parallel to the magnetic field. Such an ordered phase is the ferromagnetic phase accompanied by long- 
ranged order. Namely, the spatial correlation of two spins does not vanish when their distance approaches infinity. 
If the temperature is sufficiently high, spins strongly fluctuate in space and the average magnetization is zero. Such 
a disordered phase is the paramagnetic phase. In this phase, a nonvanishing average magnetization can be triggered 
only by an external magnetic field. Unlike the ferromagnetic phase, the spin-spin correlation decays exponentially 
when the distance exceeds the correlation length. At some critical temperature a transition from the ferromagnetic 
to paramagnetic phase occurs, accompanied by large spin fluctuations. 

This scenario has an analogy in the supersymmetric nonlinear a model. As mentioned above, the supermatrix Q 
field mimics the spin. In high dimensions {d > 2), the metallic phase - the 'ferromagnetic phase' - is formed, if 
disorders are weak enough (large vDq), in which the 'spins' Q tend to align. In the presence of Cj they all point 
towards 'north pole', A. 'Spin' fluctuations in space characterize diffusive modes (diffuson and cooperon) and their 
interactions, and long-range order is manifested in the spatial extension of wave functions. In this ordered phase, 'spin' 
fluctuations are weak characterizing a metal suffering weak localization. If disorder is strong enough (small coupling 
constant go) or in low dimensions (d < 2), the system is in the insulator phase - the 'paramagnetic phase' - where 
the 'spin' Q fluctuates strongly in space. In this disordered phase, the correlation of the 'spin' Q falls exponentially 
at large distances with the 'correlation length' being the localization length ^. This is manifested in the exponential 



localization of wave functions. As we have seen from the scaling theory (99), in dimension d > 2, a metal-insulator 
transition is triggered by lowering the coupling constant g^. this is the ana ogy of (thermodynamic) ferromagnetic- 
paramagnetic transition. We wish to point out that the above supermatrix - spin analogy should not be carried too 
far: notably, the supermatrix Q field does not play the role of the order parameter in the Anderson transition j88| . 
In fact, the average 'magnetization', /cAQ, does not vanish in both phases. Physically, this is so because the density 
of states is not critical in the Anderson transition. 

(For simplicity, below we focus on the static case, i.e., w = 0.) Open media (without internal reflection) differ from 
infinite media in that outside the medium, the 'spin' Q is fixed, pointing towards 'north pole' A (the red arrows in 



Fig. 10 1. Due to the 'ferromagnetic' coupling, near the boundary the 'spin' Q tends to align with A. As the distance to 



the interface increases, the Q-field fluctuations becomes stronger (Fig. 10 1. Since such fluctuations characterize wave 
interference, the interference strength increases as waves propagate from the interface to the sample center. Since 
homogeneous wave interference effects lead to a spatially homogeneous diffusion coefficient, D[lj), in infinite media, 
one may naturally expect that in open media a local diffusion coefficient D{r]Co) results. 
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2. Local diffusion equation and single parameter scaling of D{r;uj) 



Let us first consider the boundary constraint. The role of the second term of Eq. (78 1 is to aUgn Q to A. Indeed, if 
we ignore for the moment the first term and substitute into it the parametrization (83), we find W{y G C) = 0, i.e., 
Q{y € C) = a. For C, much smaller than the inverse of the typical value of the boundary derivative QV„(r)Q which 
is ~ C'(min(i, ^)), we may keep the leading order expansion in obtaining 



(CV„(r) - l)W = 0. 



(101) 



This implies that - if extended to the air from the random medium - would decay exponentially for distances to the 
interface much larger than C- In other words, the Goldenstone modes tend to penetrate into the air for a distance ~ C,. 
This translates the canonical physical interpretation of the extrapolation length into the field theoretical language. 
More precisely, for C <C min(^,L), the boundary constraint (78 1 is simplified to 



Qlrec = A, 



(102) 



where C is the trapping plane outside the medium and of a distance C to the interface C. Then, effect of the simplified 



boundary condition ( 102 1 is to enlarge the size of the medium by a thickness of C, which is well known for diffusive 
media |104j . 



It should be emphasized that the boundary condition ( 102 ) also applies to the localized system, provided the 
condition C ^ ^ is met. This implies that the common wisdom 104^ for a diffusive (slab) sample ~ that its effective 
thickness is larger than its actual thickness by - can be extended even to a localized sample. Most importantly, 
such common wisdom should not be carried to the case of strong internal reflection where C, <^ min(^, L). In this case. 



the perturbative expansion for Eq. ( 78 ) ceases to work and Eq. ( 102 1 breaks down. 



In this review, we are interested in the simplest case in which internal reflection is absent. In this case, C, ^ I is, 
much smaller than any other macroscopic scale. For this reason, we may ignore the slight difference of 2C, between the 
effective and actual thickness of the sample, and from now on set C — C. W e then repeat the procedures of Sec. V A 
for deriving Eq. ^. A summary of the essential differences from Sec. |V A| in the final results is as follows: (i) the 
propagator 3^o is no longer translationally invariant, i.e., 3^o(r, r'; w) ^ 3^o(r — r',0;w). The translational symmetry 
is broken by the boundary condition: 3^o(r G C, r';a)) — which is a result of Eq. (102). (ii) Because of (i), the 



probability density 3^o(r,r;w) is no longer homogeneous, but rather, depends on the position r explicitly, (iii) Up to 
two loops the weak localization correction is 



<5i?(r;w) =i?o 



3^o(r,r;w) 



C 



(103) 



which as a result of (ii), is no longer homogeneous. Here, C is zero (unity) for GOE (GUE) systems. The first term 
was also obtained in Ref. [5S by the refined diagrammatic perturbation theory, (iv) Remarkably, the correlator y 
solves Eq. ^ instead of Eq. ( 2| , with the local diffusion coefficient given by 



D{v,Cj)^Do + 6D{v- 



(104) 



For r being at the sample center, sending L to infinity (so that effects arising from interfaces do not play any role), 
we recover ordinary single parameter scaling equation (99) from Eq. (103) for d near the critical dimension. 



Most strikingly, Eq. ( 103 ) shows that the local diffusion coefficient depends on r (as well as a)) via the scaling factor 
A(r;a)) = {■Kv)~^yQ{r, r; Co). In fact, such a dependence persists to higher orders. This shows that the local diffusion 
coefficient exhibits novel single parameter scaling, i.e.. 



Do 



(105) 



where £'oo(A) is some scaling function depending on both symmetries and dimensions. It is important that the novel 
scaling (105) is valid also for large A, the non-perturbative regime. This scaling, as we will see in the next section, is 



the key component in the description of macroscopic wave transport through one-dimensional open localized media. 
It is missed by the SCLD model [Slj. We emphasize that the inability of the SCLD model to capture this scaling is 
intrinsic to the phenomenological assumption of one-loop self-consistency in open media, but not to the VW theory 
[721 [75] which is formulated for infinite systems. Indeed, if we set r to the sample center and send L to infinity, the 
scaling factor diverges as In a) ((I>~^/^) in two (one) dimensions, and the VW theory captures these infrared divergences 
and was designed for the purpose of (partly) summing up these divergences. 
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FIG. 10: The supermatrix Q and the supersymmetric nonhnear a model mimic the spin and the Heisenberg model of ferro- 
magnetism, respectively. The 'spin' Q is frozen on the transparent interface (red), and fluctuations increase inside the medium 
(blue). 



The diffusion coefficient ( 104 1 and tlie local diffusion equation ^ have important implications. (For media without 
internal reflection,) locaUzation effects are absent on the interface, i.e., D(r £ C]Cj) = Dq, which is robust against 
wave interference. In particular, there is no boundary layer which scales in the same way as the diffusion constant in 
an infinite medium, as was assumed by the earlier phenomenological model |61| . As indicated by Eq. ( 103 1, this arises 
from the inhomogeneity of wave interference effects in open media. Near the interface, wave interference characterized 
by 3^0 (r, r; w) is largely diminished due to wave energy leakage through the interface. 



3. Origin of the novel scaling 



It should be stressed that the novel scaling (105) applies to the local diffusion coefficient and therefore describes 
the unconventional hydrodynamic behavior of wave propagation in open media. It is fundamentally different from 
the ordinary single parameter scaling theory of the Thouless conductance |45j, which provides no information on 
how wave energies flow from one point to the other inside the medium. As we will see in Sec. |VI C[ the former 
leads to the latter and contains significantly more information on wave transport. Now, let us discuss the origin of 
the novel scaling. First of all, optical paths may self-intersect, forming a loop. It is well known that because of the 
time-reversal symmetry two optical paths may counterpropagate along this loop and interfere constructively with each 



other, leading to the (one-loop) weak localization correction |105l \106\ homogeneous in space (Fig. 11 upper left). 
In the presence of open interfaces, the interference picture is modified. That is, wave energies leak out the system 
through the interfaces (Fig. [Tl| upper right) and as such, at time t the probability for a path to return to the vicinity 
(of size ^ X'^~^dt) of its departure point explicitly depends on the departure point r: the returning probability is 
^ 3;o(r, r' = v;t)X'^-^dt and increases as r moves into the medium. Here, 3^o(r,r';0 is the inverse Fourier transform 
of 3^o(r,r'; w). Since wave interference enhances the backscattering probability, the diffusion coefficient is suppressed 
by an amount 



r;a; 



X^-^ I dte'"'yo{T,v-t) ~ -A(r;w). (106) 



As waves penetrate deeper inside the medium, more complicated wave interference effects take place. For example, 
optical paths tend to return to the vicinity of its departure point more (say n) times (e.g.. Fig. |11[ lower panel) 
[fSHUBa, with a probability - (A-^-^)" d^idia • • • dtn^iiVolr, r; t - ii):yo(r, r; ti - ia) • • • 3^o(r, r; t„_i). 
Two paths taking these loops as their routes pass the n loops in different orders. Having (almost) the same phases, 



33 



they also interfere constructively. This results in a higher order wave interference correction, 
5£>(")(r;£l;) 



Do 



{X'^-'r I dte^'^' [ dhdt2---dtn-iyoir,r;t-ti)yoir,v;h-t2)---yoir,r;t„_i) 

J Jti>t2>--->t„_i 

{X{r;Co)r. (107) 



Thus, the position r enters into all the wave interference corrections via the position-dependent return probability, 
which justifies the novel scaling. It is important that for n > 2, the two optical paths, though passing the loops in 
different orders, may propagate along every loop in the same direction (clockwise or anticlockwise). Such constructive 
interference picture does not require time-reversal symmetry, and is key to establishing both the local diffusion and 
novel scaling in GUE systems. 




FIG. 11: Left: examples of interference picture in infinite media. Right: in open media, wave interference is suppressed because 
wave energies leak out of the system (dashed line). Upper: interference picture underlying one-loop weak localization. Lower: 
examples of more complicated interference picture. The green line, C, stands for the air-medium interface. 



4-. Unconventional Ohm's law 



The existence of the local diffusion coefficient leads to an interesting phenomenon - the unconventional Ohm's law 
in the steady state (w = 0). Indeed, if the transverse plane of the slab is infinite, the local diffusion coefficient is 
homogeneous in the transverse direction. As a result, D{r;Cj = 0) = D{x;uj — 0) = D{x), where x is the longitudinal 
coordinate x. In the steady states, the system may sustain a uniform macroscopic current in the longitudinal direction, 
j, going through the slab. An intensity profile across the sample, I{x), which is uniform in the transverse direction is 
thereby built up. Eq. Q gives the local Fick's law, 

J = -D{x)d^I{x), (108) 

according to which the intensity bias across the sample is I{L) — /(O) = ~j Jq TJ^- Since the ensemble-averaged 
transmission {T{L)) oc j/{I{L) — /(O)), we find 

inL))-'- (109) 
Jo D{x) 

It suggests that the total resistance of the system, {T{L))^^, is the integration of the infinitesimal resistance, j^^- 
In this sense, it is similar to the familiar Ohm's law. However, it should be emphasized that in sharp contrast to the 
ordinary Ohm's law, the resistivity D^^{x) is no longer inhomogeneous in space. In particular, we will see in the next 
section that for localized samples D^^{x) may increase by many orders of magnitude as the position changes from the 
interface to the midpoint, and it is the dramatic enhancement of D~^{x) near the sample center that leads the system 
to exhibit a global localization behavior namely the exponential decay of the average transmission in the sample 
length. (For diffusive samples, such position-dependence is weak and does not lead to any interesting phenomena 



other than the scaling {T{L)) ~ 1/i.) For this reason, we term Eq. (1091 'unconventional Ohm's law'. Such law was 
first conjectured in Ref. |64) . The microscopic justification of Eq. Q puts this important conjecture on a firm level. 
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VI. LOCAL DIFFUSION IN ONE DIMENSION 



In the one-dimensional case, the asymptotic behavior of the scahng function I?oo(A) at large A has been found 
analytically [6^, which will be discussed in this section. In this case, the macroscopic diffusion equation Q is 
simplified to 



{—iCj — dxD{x; uj)dx)y{x, x'\Cj) — S{x — x'), 
which is implemented with the boundary condition, 

y{x,x']ib)\x=OorL = 



due to transparent interfaces. Single parameter scaling ( 105 1 is simplified to 

D{x;u)) 



Do 



D^{\{x;Lj)), 



and the scaling factor to 



A(x; ui) 



yo{x,x;uj) 



L coshvtJ* — cosh(vw*(l — 2x/L)) 



(110) 
(111) 

(112) 
(113) 



with uj* = —ULL'^/Dq. Throughout this section the 'effective' local density of states is the product of the local density 
of states and the cross-sectional area of QID samples, S, and we denote it as v also in order to simplify notations. 
The localization length ^ ^ vOq [91]. We emphasize that Eqs. (110l-(113) are valid for both GOE and GUE systems. 

Below, we will study in details the static case for which lu = 0. We will present analytic results for the local diffusion 
coefficient, D(x) = D{x,uj — 0) and then compare them with the results obtained from numerical simulations. As we 
will see, the analytic prediction of D{x) is entirely confirmed by numerical simulations. Thus, the novel scaling of the 
local diffusion coefBcient fully captures disorder-induced high transmission states (resonances) dominating long-time 
wave transport. As such, the seemingly contradictory properties of wave propagation, namely the local character of 
macroscopic diffusive behavior and non-local character of resonances, are unified. Here, by 'local' it means that Pick's 
law is valid everywhere in space, while by 'non-local' that resonances result from that waves propagate back and forth 
between two interfaces. 



Static local diffusion coefficient 



For Cj — 0, Eqs. (112) and (113) arc further simplified to 

D{x)/Do ^ D^{\{x)), \{x-uj - 0) = \{x) = x{L - x)/{L^). 
To proceed further, we introduce a new scaling function, 

g{\) = \-^DU\). 



The single parameter scaling (114) gives the Gell-Mann-Low equation, 

d\ng 



dlnA 



with the /?- function depending only on Q. 

For ^ 1, according to Eq. (1031 the /3- function is expressed in terms of the l/t/-expansion. 



/3(^) = -l 



Cl 

g 



C2_ 



e > 1. 



(114) 



(115) 



(116) 



(117) 



For GOE (GUE) systems, the nonvanishing subleading term is ^ g~^ {g~^)^ and the corresponding coefficient ci (02) 
is negative which arises from the inhomogeneous leading weak localization correction. As such, Eqs. (116) and ( |117 ) 
bear an analogy to the ordinary single parameter scaling theory ([99]) : the scaling function g here plays the role of the 
'Thouless conductance', 5, and A of the 'size', Lq. Identifying this analogy, one may follow Ref. flB^ to extrapolate 
the Gell-Mann-Low function into the regime of small g and obtain 



p{g) = \ng, g«i. 



(118) 
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In fact, this non-perturbative /3-function has been justified by using the exact solution of the correlator ( 81 1 for special 
coset space leading to Q € GMat(3, 2|A) [95 t I159j . On the other hand, it has been well established ^88^^95* iSZj that 
provided strong QID localization exists, its properties are insensitive to the symmetry of Q. The latter only affects the 
numerical coefficient in the analytic formula of the localization length j941 195j . Since introducing this exact solution 
requires very advanced mathematical knowledge namely Fourier analysis on a hyperbolic supermanifold with constant 
curvature [95H97j which is far beyond the reach of the present review, we shall not discuss it further. Instead, we 



will present the numerical confirmation of the non-perturbative /3- function (1181 as well as a straightforward physical 
interpretation below. 



Notice that single parameter scaling (114 1 is valid for all the sample lengths (larger than the transport mean free 
path). Therefore, we may send L — >■ oo and the sample becomes semi-infinite. In this case, the scaling factor A has a 
simple form of x/^. On the other hand, the probability for waves to to penetrate into the sample of a distance a; 3> ^ 
is exponentially small ^ e^-^/^ [48 ST. Therefore, D{x) decays exponenti ally f or x ^ ^, i.e., D{x)/Dq ^ e~^l^ . This 



leads to a large- A asymptotic scaling function, D^iX) ~ e . Equation (118) is thereby reproduced and thanks to 



the universality of Eq. (114), is independent of the sample length. 



The novel single parameter scaling theory namely Eqs. (116) and (118) - for the local diffusion coefficient - gives 
a very remarkable prediction for the local diffusion coefficient deep inside the localized sample, namely Eq. ([6|. The 
result (solid lines in Fig. |4| is contrary to the result of SCLD [64]: instead of simply decaying exponentially from the 
interface (dashed lines in Fig.|4|, inside the sample D{x) acquires an enhancement factor of e^^/(^^\ which drastically 
reduces the falloff of D(x) into the sample. Such a dramatic enhancement is also valid in one dimension where the 
localization length is of the order of the transport mean free path. 

B. Numerical evidence of the novel scaling 

Numerical experiments on the spatially resolved wave intensity across a randomly layered medium were performed 
in Ref. |69| by adopting the method of Ref. |21j . The random medium with layer thickness a is embedded in an 
air background. The relative permittivity at each layer fluctuates independently, and is uniformly distributed in the 
interval [1 — cr, 1 + cr] (the so-called 'rectangular distribution'), where a measures the degree of randomness of the 
system and is set to 0.7 in numerical experiments. The air-medium interfaces are transparent. Then, a plane wave 
of (angular) frequency w is excited, and two frequencies, w = 1.65c/a and 0.72c/a, are considered, where we have 
restored wave velocity in air. For each w, two million dielectric disordered configurations are realized for different 
sample lengths. Then, for each disordered configuration, the standard transfer matrix method is used to calculate 
the transmission coefficient T and wave field E(x), from which one obtains the ensemble-averaged current j = (T) 
and wave intensity profile I{x) = (|i?(x)p). Since the current across the sample is uniform in steady state, by further 



presuming the local Fick's law namely Eq. (108) one may compute D(x). 

The numerical results of the local diffusion coefhcient are shown in Fig. |4] (squares and circles) for five different 
sample lengths. In the limit L — > oo, D{x) decays exponentially from the interface (dotted line). This confirms 



the large ^/-asymptotic of the /3-function namely Eq. (118). Most importantly, data points for different frequencies. 



1.65c/a (squares) and 0.72c/a (circles), overlap which indicates that the scaling behavior is universal independent of 
the microscopic parameters of media. The universal curves are in excellent agreement with those predicted by Eq. ([6]) 
(solid lines). Obviously, except in the small regime near the interfaces, the results from simulations are significantly 
larger than those obtained by the prediction of the SCLD model (dashed lines) and the deviation is more and more 
prominent as the ratio L/£^ increases. 

C. Roles of high transmission states 
Dynamic studies of wave transport through localized samples [H] have shown that the long-time transport in these 



systems is dominated by rare long-lived modes with high transmission coefficients (see Sec. IIBl). This indicates 
that in the strongly localized regime far-reaching consequences may arise from the Lyapunov exponent fluctuations. 
To see this let us start from a semi-quantitative analysis. A basic prediction of the SCLD model is as follows that 
for localized samples D[x) is a sum of two symmetric exponential decays truncated at the sample mid-point (see the 
dotted line in Fig. [i] which represents the left truncated exponential decay) . This result may be written as 

/•CO 

D{x)^ d7i^(7)e-^'-"(-'^-^), P(7)=<5(7-r') (119) 



for X sufficiently away from the interfaces, i.e., min(x, L — x) ^ ^. Importantly, the above distribution P(7) reflects 
the absence of fluctuations in the Lyapunov exponent. But this is correct only if the sample length is infinite. (That 
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FIG. 12: Local diffusion breaks down upon removing a small portion of high-transmission (T > Tc) states from the original 
ensemble. The source is placed at a; = 0. (from Ref. [BS] with reproduction permission from C. S. Tian, S. K. Cheung, and Z. 
Q. Zhang © The American Physical Society) 



is, a semi-infinite sample is considered.) For finite-sized localized samples, the distribution P(7) is given by Eq. ([S]) 
instead. Replacing the Dirac distribution in Eq. (1191 with the distribution ([8| we reproduce the correct result for 
D{x) namely Eq. i.e., 



Dix) 







-L/(45) 



(120) 



This clearly shows that the SCLD model completely ignores the fluctuations of the Lyapunov exponent. Interestingly, 
the second line in Eq. (120) shows that the integral is dominated by 7 '--^ |1 — 2x/L\£^~^. In particular, for x closed 
to the sample center, the integral is dominated by the small 7 tail of the distribution -P(7). This suggests that the 
significant enhancement of the local diffusion coefficient from the exponential decay near the sample center is deeply 
rooted in the rare high transmission states. 

In Ref. [SS] , the roles of high transmission states in the static transport are studied numerically, and data leading 
to results shown in Fig. |4] are further analyzed. First of all, InT follows the normal distribution ([s]) [IHl [SO]j with 
the average ~ —L/^ and variance ~ 2L/^. Then, to investigate the role of high transmission states a small portion 
of them (with transmission T > Tc) are removed from the original ensemble (w = 1.65c/a and L/^ = 10). Most of 
removed states are singly localized states, and the remainder (with a small portion) are necklace states |164j . D{x)/Dq 
is re-computed for the new ensemble. As shown in Fig. [12] even when the fraction of states removed is as small as 
0.6% (solid line, in red), the result deviates drastically from the original one (solid line, in black) and is asymmetric. 
This signals the breakdown of local diffusion and indicates that rare high-transmission states are intrinsic to local 
diffusion and novel scaling. 

Substituting Eq. ^ into Eq. ( |109[ ), we find 



{T{L)) 



(121) 



We see that the dramatic enhancement of the resistivity near the sample center leads the system to exhibit a global 
localization behavior namely the exponential decay of the average transmission in the sample length, as we mentioned 
above. On the other hand, from the normal distribution of ^(7) it follows that 



d{\nT{L)) ^ _ 1 
dL ^ ' 



(122) 



Eqs. (|12lj) and |122[ ) show that (for both GOE and GUE systems 
and geometrica 



means of T{L) differ by a factor of 4 
transmission, and previously was discovered for the conductance distribution 



the localization lengths obtained by the arithmetic 
This factor arises from the rare disorder-induced resonant 
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VII. CONCLUSION 



How do classical electromagnetic waves propagate through random open media? This is a long-standing fundamental 
issue in Anderson localization. One might expect that like many complex systems propagation of waves in random 
media bears the microscopic and macroscopic description simultaneously. The former is built upon the eigenmodes of 
the Maxwell equation which provide complete information on wave scattering (hence the term 'microscopic'), while 
the latter describes the system in terms of a single macroscopic variable - the average energy density. It satisfies some 
macroscopic (generalized) diffusion equation which is valid only on the scale much larger than the mean free path (hence 
the term 'macroscopic'). For open media, the microscopic description of waves can be traded to the superposition of 
excited quasi-normal modes. In contrast, the rational of the macroscopic description especially for strongly localized 
open media has been an open question explored by many researchers and had largely not been established as yet. This 
review is devoted to substantial recent progress achieved in the macroscopic (or 'hydrodynamic' called by condensed 
matter physicists) description of wave propagation in random open media. 

We reviewed the recently developed first-principles theory for classical wave localization in open media |67L I69j that 
justifies unconventional macroscopic wave diffusion in these systems. In essence, it describes localization physics in 
terms of an effective action of the supermatrix Q field. The theory differs crucially from the supersymmetric field 
theory of (electron) localization in infinite systems in the air-medium coupling action introduced by wave energy 
leakage through the interface. The latter breaks the translational symmetry of the low-energy field theory and 
constrains the supermatrix field by some boundary conditions. We showed at a pedagogical level how to use this 
theory to explicitly calculate the wave intensity profile and to establish a macroscopic description of wave propagation 
at a firm level. 

We reviewed highly unconventional macroscopic wave diffusion in open media. This is described by a linear macro- 
scopic equation - the local diffusion equation. It differs from the normal diffusion equation in that the diffusion 
coefficient is inhomogeneous in space. Most importantly, it was discovered in Ref. [52] that the local diffusion coeffi- 
cient exhibits novel single parameter scaling. That is, it depends on the distance to the interface via a scaling factor 
proportional to the returning probability density, and the universal scaling function depends on both the system's 
symmetry and dimensionality. In the static case (steady state), the scaling function in one dimensions has been found 
analytically giving the explicit analytic expression for the static local diffusion coefficient. These results are fully con- 
firmed by numerical simulations. This suggests a profound new concept which is contrary to previous expectations. 
That is, the resonant transmission - formed in the scale of system's size - does not wash out Pick's law (thereby 
macroscopic diffusion) which is valid locally in space; rather, it results in novel scaling of the local diffusion coefficient 
which plays a decisive role in the long-time transport of localized waves. 

From the practical view point, the novel scaling of the diffusion coefficient (1051 and the local diffusion equation 
Q are sufficient for theoretical analysis of many experimental measurements involving the averaged wave energy 
density. Such approach based on the macroscopic equation has the great advantage of technical simplicity over its 
sophisticated microscopic theory. In this case, one has to find the scaling function of the local diffusion coefficient by 
other methods (e.g., numerical simulations or fitting experimental data), and this is a relatively simple task. 

We wish to mention some important directions in future studies. While the novel scaling of local diffusion coefficient 
has been shown to generally exist in the slab geometry in arbitrary dimensions, the non-perturbative scaling function 
remains to be worked out explicitly, and this could have far-reaching consequences. In fact, the dynamic local 
diffusion coefficient is related to a number of important issues. These include time-resolved transmissions [191 |2T], the 
statistics of quasi- normal modes [31] , dynamic single parameter scaling [21] , anomalously localized states and dynamic 
conductance [140j, etc.. Furthermore, the local diffusion coefficient (both dynamic and static) in high dimensions may 
help to reveal new physics due to the strong interplay between the openness of the medium and criticality. Indeed, we 
have seen that the novel scaling already leads to highly non-trivial results in one dimensions. That is, the localization 
length given by the mean transmission is four times larger than the typical localization length. We naturally expect 
that its effects in high dimensions might be even more interesting. Whether and to what extent the novel scaling of 
local diffusion coefficient affects the size-dependence of the conductance near Anderson transition? How does it lead 
to the pulsed wave decay on the output interface? These questions are fully open. 

In this review, we have focused on ideal optical systems ignoring internal reflections at the boundary and absorption 
(gain) in the bulk. In experiments, these can have significant effects. Interestingly, there have been experimental 
and analytic evidence |166j showing that the common wisdom in optics - adding the extrapolation length to the 
sample length |104j - breaks down. The interplay between internal reflection and wave interference may lead to 
very rich localization behavior. The other issue, effects of (linear) absorption and gain in random media has received 
considerable attention (51 [71 [551 11671 1168] . but these works deal with inflnite media. As far as realistic optical devices or 
experimental environments are concerned, one often deals with open media and this will lead to intriguing phenomena 
(for examples, see Refs. [T8 l 11411 1142[ 11691 1170 |). A fundamental problem therefore is how these factors interact with 
localization and unconventional macroscopic diffusion. 
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The supersymmetric field theory reviewed here is formulated for scalar wave systems. As mentioned above, more 
complete studies of localization of classical electromagnetic waves will require a formulation for vector wave systems. 
What will be the symmetry of the supermatrix Q and its low-energy action then? To what extent will the vector 
nature modify the macroscopics of localized waves in open media? These important questions are central to studies 
of light localization and largely unexplored. Throughout this review, we have focused on classical wave systems. In 
principle, the macroscopics reviewed here exists in electronic systems as well. As far as the latter is concerned, it 
is well known that the Anderson transition is classified into ten symmetry classes |171j . Of particular interests is 
the sympletic class. For this symmetry class, an intriguing phenomenon intrinsic to the openness of the system was 
discovered some time ago [HZ]. That is, even though an infinite QID disordered wire exhibits strong localization |94j . 
the localization behavior is dramatically changed as long as the system is open (namely coupled to ideal leads): the 
conductance decreases to a constant instead of zero as the sample length goes to infinity. It is therefore conceivable 
that the symmetry may affect profoundly the behavior of the local diffusion coefficient. 

Finally, we emphasize that throughout this review we consider only ensemble averaged observables. The theory 
reviewed here provides no information on individual disorder configurations. The resonant properties of individual 
disorder configurations in the localized regime can be studied by other approaches (see, e.g., Refs. [5TII1721 1173| and 
references therein). 
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Appendix A: The basics of Grassmann algebra and supermathematics 



This Appendix includes some preliminary definitions and basic theorems of Grassmann algebra and super- 
mathematics which are required for following (most) technical details of this review. We do not aim at a 
complete introduction to mathematical foundation of the supersymmetric field theory, for which we refer readers to 
Refs. [88, 93, 139, l'74Hl'79l. 

Grassmann algebra. The Grassmann algebra A (over complex C) is an algebra constructed from n generators 
Xi- These generators - the so-called Grassmannians - satisfy the anticommuting relation, i.e., 

XtX] + X]Xi = 0, Vi, j. (Al) 

From this definition it follows that = if the integer k > 2. This gives an important theorem: all elements 
in the Grassmann algebra A are first degree polynomials in each generator Xi. In other words, a generic function, 
/(xi,X2, • • • ,Xn)i must be constructed from 2" monomials, Xi^X2^ ' ' 'Xn") where the power ak takes the value of 
either or 1, i.e., 

f{Xl,X2,--- ,Xn)=Y^ Ca,a2-a„Xl'xT ---Xn", (A2) 

with the coefficient Ca^a2---a„ G C. 
Define the reflection: 



P{X^) = -X^■ (A3) 

Then, a monomial of degree k has a parity (— I)'' for P{xiiXi-2 " ' Xik) = (~l)'^XiiXi2 ' ' ' Xi*, ■ This divides the 
Grassmann algebra into two parts, A = A'^ U A^ , A'^ n A^ — 0, where A'^ (A^) has even (odd) parity, constructed 
from monomials of even (odd) number degree and therefore composed of commuting (anticommuting) elements. The 
former is a subalgebra in A whose elements are similar to ordinary complex numbers. 
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The complex conjugate operation (adjoint) of the Grassmann algebra is defined by 

ixi)* = xh iXiT = -Xi, {xiXjT = xlx*j- 



(A4) 



The '— ' sign in the second relation makes xtXi act like a 'real number', i.e., {x*Xi)* = xtXi- Note that x* is 
independent of Xi and as such, the complex conjugate of the Grassmann algebra is purely a formal definition. 



Calculus in the Grassmann algebra. Eq. (A2| allows us to introduce the partial derivative with respect to 



the generator Xi- Because any elements in the Grassmann algebra A, which are functions of the generators, can be 
written as Eq. ( A2 ) , we may move Xi to the left for each term, obtaining 



f = fo + X^fl, (A5) 
where the coefficients /o,i are Xi-independent. Then, the left partial derivative with respect to Xi is defined as 



dXi 



h 



(A6) 



(The right partial derivative can be introduced in the similar way.) It is important that in moving Xi to the left, a 
fermion sign is left due to the anticommuting relation ( Al ). It is easy to check that the partial derivative defined by 
Eq. (A6) is a nilpotcnt operator, i.e., {d/dxif' — 0. Furthermore, if we consider Xi as an operator acting on A via 
left-multiplication, then the nilpotcnt operator d/dxi and the operator Xi constitute the Clifford algebra. 



XiXj + XjXt = 0, 



d d 9 9 „ d d ^ 

= 0, X^^ 1-^— Xi=0»J- 

dXj dxj 



dXi dXj dxj dx, 



(A7) 



For Grassmann algebra, the integration operation is defined as the (left) derivative. More precisely, the single 
variable integration is defined as 



while the multiple variable integration as 



fdXi 



fdxk ■■■dxi 



d 
dxi 



d 



d 



dxk dx) 



(A8) 



(A9) 



In fact, this formal definiti on sa tisfies the properties of the ordinary definite integral. To see this let us denote the 
integration defined in Eq. (A8) as I. First of all, it is obviously a linear operator, I[fi + = I[fi] + I[f2] for 
/i,2 S A. Secondly, because the partial derivative is a nilpotcnt operator, we have I[df /dxi] = implying that the 
integral of a total derivative vanishes and dl[f]/dxi = implying that a definite integral is a constant. Thirdly, if 
/' is Xi-independent, we have /[//'] = /[/]/': the overall factor can be pulled out of the integral (from the right). 
Finally, if we make the replacement Xi ^ Xi + 'yij with r]i a constant Grassmannian, the integral remains the same 
implying dxi = d{xi 

Gaussian integrals with Grassmannians. In this review we are interested in the case where the generator Xi 
and its adjoint x*i (* = 1? 2, ■ ■ • , n) appear in pairs. Consider an invertible nx n matrix, M, with the entries Mij G C. 
By Taylor expanding e~-^>*^'3Xj ^^d using Eq. (A9), one finds an important identity. 



J e-^'^'^'^^dxldxi ■ ■ ■ dxldxn = detM. 



(AlO) 



This is in sharp contrast to the Gaussian integral over ordinary complex variables, S* , S'i, i = 1, 2, • • • , rt, which is 

-S'MijS, dS^dSi ^ ^ ^ dS*^dS.. 



(detM)" 



(All) 



Importantly, the difference of the power of the determinant, ±1, on the right-hand side of Eqs. ( |A10 ) and (All ) refiect 
different rules of the change of variables for Grassmiannians and complex variables. Indeed, the linear transformations: 
X' = Mx and S" = MS lead to dx« • • • dxi = det M dx'n ■ ■ ■ dx'i while dSn--- dSi = (det M)'^ dS'^--- dS[. By using 
Eq. (AlO I, we obtain another important identity: 



(detM)-^ / x^XJe 



~X:M.,X,dxldxi ■ ■ ■ dXndXn = {M-% . 



(A12) 
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The supervector and the supermatrix. A supervector </> is composed of an n-anticommuting component vector 
X and an m-commuting component vector S, (We shall consider the case of n = m.) 



s 









f 


X = 


[Xnj 


, s = 


[sj 



(A13) 



where Xi i^i) s-re anticommuting (commuting) variables. All such supervectors constitute a linear space, The 
transpose of the supervector (/> is defined by 



iXl ■ --XnSl 

while its Hermitian conjugate by 

<P^^{xl---XnSV 

The inner product on the supervector space 3>n,n is defined by 



(AM) 
(A15) 



') 



(A16) 



The norm -y/ (0, (f)) defines the 'length' of a supervector </>. 

The linear transformation on the supervector space is described by a supermatrix, M, which transfers anti- 
commuting (fermionic) into commuting (bosonic) components and vice versa. It has the following structure, 



M . 



Mff Mfb 
Mbf Mbb 



fb 



(A17) 



with the superscript 'fb' standing for the fermion-boson space accounting for the supervector structure defined in 
the first relation of Eq. ( |A13 1. Here, A/ff, AfeB, M^b^ and Mbf are n x n matrices. Furthermore, the entries of 
diagonal (off-diagonal) blocks, A/bb, A^ff (A^bf, A/fb), are commuting (anticommuting) variables belonging to ^+ 
{A~). The supermatrix product is defined as usual. 

The transpose of the supermatrix is defined by {M(j))'^(p' = (fp- M"^ (f/ , V0, (f)' G $„.„. This gives 



AfT = 



A//Tp 



-<F 



Mfb ^bb 



fb 



(A18) 



where the '— ' sign in the upper right block arises from the anticommuting relation, and for each block Maa' the 
transpose is defined as usual. By the definition one may easily check that for two supermatrices, Af, Af', 



{MM 



M'^^M"^ 



(A19) 



{M(t))^ = (f^M'^ 



(A20) 



The Hermitian conjugate of the supermatrix is defined by (M0, 0') = (0, Af^0'), V0, 0' € $n,n- This gives 
M'^ = (Af""")*. The pseudounitary transformation can be defined as such U that conserves the me tric t ensor 
i.e., {U(t),KU(j)') EE {(j),K(j)') or WKU = K. If if = 1, then U becomes unitary. Using the identity ( |A19[ ) and the 
third relation in Eq. ( A4 ) , we have 



{MM')"^ = Af'tA/t. 



(A21) 



By using Eq. (A18) and the second relation in Eq. (|A4|, we have 

(A/^)^ = Af, 



(A22) 



which should be contrasted with (Af ) = (T3 Af crj 



fb 
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For the supermatrix, we introduce the supertrace, denoted as 'str', 

strM = trAfpF - tr7\fBB, (A23) 

where 'tr' stands for the normal trace. (This differs from that defined in some other references, e.g., Refs. [5^ I138j . 
by a minus sign.) By using the definition, it is easy to check the following relations, 

str(MM') = str(M'M) (A24) 

and 

(0, M(j)') = -str(M0' ® 0^). (A25) 

The superdeterminant is defined by 

In sdetM = str In M, (A26) 
which is the direct generalization of the well-known identity Indet(-) — trln(-). This gives 

sdetM = det(AfFF - AfFBMBB^-^BF)(detAfBB)"^ (A27) 

By the definition we have 

sdet(A//Af') = sdet AfsdetAf'. (A28) 



By using Eqs. (A18l and (A27), we can show 

sdet M'^ = sdet M. (A29) 



Gaussian integrals with supervectors. Given the supermatrix M, we have the following Gaussian integral over 
the supervector. 



= sdetAf, 



(A30) 



where the measure 



nn 7*7 TT"i dS*dSi i 
i=i dXt U^=l and 



(sdetAf)-i 



(Af"i) 



(A31) 



In Eq. ( A31 1, if the pre-exponential factor is a product of 0's with the number > 2, then the ordinary Wick theorem 
holds. The only difference is that exchanging Grassmannians leads to a fermionic sign. It is very important that here, 
we have assumed the convergence of the integral, i.e., ReS'^AfeB'S' > 0. 



:---x----x--x--» 
^ 



C 









( 

. ^ 


----X---X---X----; 


(---X---X-X--X--- 



FIG. 13: Contour deformation for the integral over the matrix element 6it,it (left) and 62*, 2t (right). 
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Appendix B: The diagonal mean- field saddle points 



In this Appendix we will find the diagonal solutions to the saddle point equation (49 1, i.e., Qo = diag{gi}. (Recall 
that the index i — mat.) Substituting Qo into Eq. (49) gives 



12 = 



1 



dkiy{\k\)- 



1 

2ioj^q, 2 
1 



'A 



dkv{k)- 



1 



(Bl) 



Notice that depending on the dimension this integral may suffer ultraviolet divergence. In this case, since the 
wavelength of photons has no lower bound, one may cure this divergence by multiplying the integrand by a factor 
g±«*C^ where the '+'('"') sign corresponds to Reg^ < (> 0), and send the positive infinitesimal 5 to zero in the end. 
(A more physical way of curing this divergence is to introduce more realistic disorders as discussed in Sec. IV A 4 In 
any case, this divergence only enters into the renormalization of the background value of the refractive index which 
is unimportant for present discussions.) In the limit Au)u{uj) <C 1 this transcendental equation is solvable, giving 



±|Regi| + ilm^i, \Reqi\ 



-Aujv{ijj) + o[Auji'{uj)), Imgi = o(Au!i'{uj)). 



(B2) 



That is. 



qi = -Aujv{Lo)\i ^ qXi, Aj = ±1. 



(B3) 



At the first glance, Eq. (B3) gives 2^ = 256 solutions. But this is not true. First of all, due to the charge conjugation 
symmetry, Q = CKQ^^ KC ^ that inherits from the reality condition (32), 



Amal — A^^ck2 



(B4) 



Secondly, for the functional integral (47) to converge, it is necessary that QmBt.mBt is purely imaginary. Therefore, 



the saddle points (B3) do not belong to this manifold. To solve this problem, in integrating out the matrix elements 
QmBt,mBt, wc dcform the integral contour such that (i) no singularities of the action F[Q] are encountered, and (ii) 
the new contour passes qiet — +q {XiBt = +1) and q2Bt = ~q (A2Bt — —1), where the '+' ('— ') sign accounts for 
the analytic structure of the advanced (retarded) Green function. To ensure the condition (i) we should first look for 
the singularities of the integrand of the partition function ( 47 ) . According to the definitions ( A26 ) and ( A27 ) , they 
are the zeros of det(— — + ww+CTf ® 1*'' + 2iu}'^QBB)- Defining Qbb — ib, we find that for the integral over the 
diagonal entries, bu^it G M, the singularities lie at the line Im6it^it = S > 0. (Recall that w"*" = ui + i6.) Taking this 
into account, we deform the contour, C, into Ci that passes through 6it,it = —iq, corresponding to Aiet = 1 (Fig. 13 
left). Similarly, for the integral over the diagonal entries, &2t,2t, the singularities lie at the line Im62t,2t = S < 
Therefore, we deform the contour into C2 that passes through &2t,2t — iq, corresponding to A2Bt = — 1 (Fig. 13 right). 
For the other cases, i.e., XiBt = ~1 and (or) A2Bt = 1, to perform the corresponding contour deformation must cross 
the singularities. Such saddle points are unphysical and must be discarded. Thus, we have 



AlBt — — A2Bt — 1- 



(B5) 



Taking Eqs. (B4) and (B5) into account, now we have only four possibilities namely Aipi — \if2 — il and 
A2F1 = A2F2 = ±1- Integrating out fluctuations around the saddle points with Aipt — X2Ft leads to a negligible small 
factor because the bosonic and fermionic degrees of freedom do not compensate. Therefore, the two saddle points 
of AiFt = A2Ft do not play any roles. Finally, only two saddle points out of 256 possibilities are left: A and —kA. 
Andreev and Altshuler discovered f l80[|181] that both are the stationary points of the action (61 ). The former (latter) 



determines the perturbative (non-perturbative) part of the level correlator in finite closed systems. 



Appendix C: Derivation of the fluctuation action 



The neighborhood of Qq — TQ^[qA)To in the supermatrix space can be parametrized by Qo + — (2o(l + 
dT))^^{qA+SA){To{l+5T)). Here 6T and 6A parametrize fiuctuations in the coset space and eigenvalues, respectively. 
Expanding 6T and SA, we obtain the fluctuation SQ ~ Tq^SATq + [Qo,i5T]. It is easy to see that the first (second) 
term (anti) commutes with Qq. Therefore, the former (latter) may be identified as JQ' ((5Q*). In other words, the 
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eigenvalue fluctuations generate the longitudinal components while the fluctuations in the coset space generate the 
transverse components. The former violates the nonlinear constraint (i.e., is a constant matrix) and thereby brings 
Qo out of the saddle point manifold. 



In the remainder of this Append ix w e derive the fluctuation (55). Upon the substitution of the decomposition 

(CI) 



<5Qq = 5Q\ + (5(3q, we simphfy Eq. (l54|) to 



rfq 

W) 



^str (H^5g^<5QLq + U\bQ\8Q'_, 



To simplify the action of the longitudinal components, we notice 



dk 

(2^ 



^0 k 



dk 



qA ^0 k-^,qA 



dk 

To 



|,T„-i(gA)r„ 



(C2) 



where in deriving the last line, we have used the fact that A = diag{Ai} with AiQ,f = ~^2at — 1- Because 



dk 

(2^ 

^ /'OO 

2 



^0 k 



q 



q ~, 



qXi] Go [k- -,q\ 



1 



dk 

{2nY 



So (k,gA.j) tJo (k, qXi) 



(C3) 



we have — A ^. Here, the same ultraviolet regularization as Eq. (Bl) was used. 
To simplify the action of the transverse components, we notice 



dk 



k 



Qo ^0 k 



Tg-Miag 



-Tg-Miag 



dk 
dk 



q ~ 



^0 k+ ^,gA, k 



(277)'* [(k + q/2)2 - a;2 + 2iuj^qX,] [(k - q/2)2 - uj^ - 2ioj'^qX,] 



Tn 



(C4) 



Since the integral over k is dominated by the shell of |k| « w, we have |k • q| ^ jqp. As a result, is simplified to 



A-1 - A-^Tn^d 
A-i/l- 



^lag' 
dQd 



)<* [k2 + k • q - w2 + 2iiu^qX,] [k2 - k • q - cj2 - 2ia;2^Ai] 
d^d 1 



1 



Tn 



n ■ q 



1 



1 



Q.d 1- 



1 



n ■ q 



d \TiAuj'^v{uj) 



dA \'kAuj'^v{ijj) 



(C5) 



where n is the unit vector in the d-dimensional unit sphere of volume ^d- 
Substituting Eq. (C5| into Eq. \C2\ gives 



Acj2;y(t^) 



dq 



(27r) 



^|q|^str{<5Q;<5QLq} 



(C6) 



44 



Upon inserting this action into the partition function (47), we may integrate out the longitudinal components, 5Ql^. 
Such a Gaussian integral is unity because the degrees or freedom of commuting and anticommuting variables com- 
pensate each other. The remaining action is a functional of the transverse components only which is 



drstr(V5Q*(r))2 



(C7) 



Since (Qo + SQ'^{r)) stays inside the manifold defined by Eq. (51), we have \/SQ*{r) 
then gives the fluctuation action (55). 



gV(r-i(r)AT(r)). Eq. (C7) 



Appendix D: Derivation of the interface action 



In thi s Appendix, we wiUderive the interface action Fi^t [Q] with the help of the Zirnbauer-Efetov theorem. Formally, 
Eq. (73) differs from Eq. (15) only in the non-Hcrmitian part of the effective Hamiltonian, i.e., zLiBSc- Therefore, we 
may repeat the procedures of Sec. |IV A"2l As a result. 



strr In I 1 



BAG, 



(Dl) 



where the supertrace 'strr' includes the integral over the spatial coordinates. Near the interface we may choose 
locally the coordinate system as r = {r±,z) where the first {d — 1) coordinates give the projection to the interface 
C and the last one the distance to the interface. Note that the z-axis points towards the random medium. As 
shown in Sec. IV A3[ (inside the random medium) the matrix Green function Qq exponentially decays in space with 
a characteristic length of the order of the mean free path. Moreover, both Q and the dielectric field vary over a scale 
much larger than the mean free path. Taking these into account, we find 



^o(r,r';Q(r')) 



TTl 



(2^) 



dk- 



cos{kz) cos(fcz') 



Kr'))2-|k^|2-P+»5^Q(r') 



r,r' e V+. 



(D2) 



Here, (cj(r))^ — w^(l + e(r)), with (1 + e(r)) the local average refractive index near the interface which generally 
causes the internal reflection. Because uj{r) locally depends on the average refractive index, from Eq. (56) we find 
that near the interface the mean free path I generally acquires a r-dependence also. 
On the other hand, the Green function g^2(r,r') is 



55(r,r') 



(27r 



dk- 



sin(fcz) sin(fcz') 



Substituting it into Eq. ( 74 ) gives 



r. r' € V_ 



/(rl) 



Let us insert Eqs. (D2| and (D4) into Eq. (Dl|). As a result 
where the coefficient 



1 .d-i 1^ (1 + (r) AQ(r)) 



^ ^2 (^l + e(r; 



(D3) 



(D4) 



(D5) 



(D6) 



In deriving Eq. (D5) we have adopted the trick used in AppendixjC] That is, we made the local gauge transformation 
so that Q{r) becomes diagonal, and after integrating out k we rotated the diagonal matrix back to Q(r). Eq. (D5) 
can be rewritten as Eq. ( 75 ) [67] , where the internal reflection coefficient is 



^kjr) 



1 - Qk^ (r) 
1 + aici (r) 



(D7) 



in agreement with the well-known result |130j . 
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Appendix E: Derivation of Eq. ( 90 ) 



First of all, we have 

F2[W] = 



T 

~2 



drstr {-2Da{ViWf{iWf + i<:j{iWf) 

drstr {2D0 {v'^iW{iWf + {ViWiWf + {ViWf [iWf^ + iCo{iWf^ , 



(El) 



where in deriving the second line we have used the integral by parts. Let us substitute it into the expansion (86). This 
leads to the leading wave interference correction to the bare correlator 53^1 (r, r'; w). The latter has five nonvanishing 
contributions, i.e., 53^i(r, r';Q) = Syia + Syu + (53^ic + <53^id + Syie- The first term is 



syia 



{strA+iW{Y)A^{iW{r')f + stYA+{iW{r)f A^iW{r'))^ 



(E2) 



corresponding to the diagram given in Fig. [s] (b) . The other terms are 



5yib = 
^yxc = 
^yid = 
?>yi. = 



4w 



4a; 



TTi/Do / dri(stryl+iM/^(r)A_iVK(r')strV2iVK(ri)(iVK(ri))3)^, 
dri {&XTA+iW{v)A^iW{x')'sX.x(iW{vx)f)^ , 



4a; 
4a; 



-kvDq J dri {atTA+iW{r)A_iW{r')stT(yiW{ri)iW{ri)f)^ , 
TTvDo [ dri{stTA+iW{r)A^iW{r')str{ViW{ri)f{iW{ri)f) 



(E3) 



given by Fig.|8](a). 

For syu let us now make the first contraction by using the rules (87) and start from the M^- factor acted by the 
operator V^. There are three possibilities, leading to 



syib = 



7ri^(a;) 
4a; 



TTvDo J drA (strA+iW{r)A^iW{r')stT\'^iW{ri){iW{ri)f'^^ 



+ {fitrA+iW{r)A^iW{r')stTW^iW{ri){iW{ri)f 



where 



sy[ 



+ (^stTA+iW{r)A_iW(r')stTV^iW{ri)iW{ri)iWiri)iW{ri)j^ \ = Sy'^ + Sy'-^^, 



(E4) 



ifc 



4a; 



Do / dri{v2(-3;o(ri,r';(I.))(str.l+iW^(r)yl_(it¥(ri))3)g 



+V^{-yo{ri,r;Co)){strA+iiW{ri)fA^iW{r'))^}, 



sy'lb 



1Il'{uj) 

4a; 



r I I 

Do j dri(stYA+iW{v)A^iW{v')siYV^iW{vi)iW{vi)iW{Yi)iW{ri)^ ^. (E5) 



Likewise, we have 



7ri^(a;) 
4a; 



mvCj / drA UtT:A+iW{Y)A^iW(T')sinW{Yi)(iW{Yi)y 



+ UiT:A+iW{Y)A^iW{Y')sinW{Yi){iW{Yi))'' 



iCj / dY^{{^yo{Yuv'-Cj)){BiYA+lW{Y)A_{iW{Yi)f) 



7ri^(a;) 
4a; 

-{-yo{Y^,Y- Co)) {sir A+(iW{Y^)fA^iW{Y'))^^], 



(E6) 
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2 



and 



sy'ld 
syie 



Auj 



TTiyDo J dri (stTA+iW{r)A^iW{r')strViWiri)iWiri)ViW{ri)iW{ri)'j ^, 

/' I 1 

dn (8tTA+iW{r)A^iW{r')stT\/iW{ri)iW{ri)ViW{Ti)iW{ri)'^^ , (E7) 



7rj/(aj) 
Auj 



2ttuDo J dri (^strA+iW{r)A^iW{r')stTViW{ri)ViW{ri)iW{ri)iW{ri) 
TTuDo J dri(strA+iW{r)A_iWir') 



,11 I 
cstr {iW{ri)V^iW{ri)iW{ri)iW{ri) + iWiri)ViW{ri)iWiri)ViW{ri 



Here, we have used the integral by parts to obtain the second equahty of Eq. (E8|. 
By using the definition of the bare propagator ( 88 1 , we find 

syia + sy'u + syic = o 



(E8) 



(E9) 



which is a resuh of the energy conservation law (cf. Eq. (93l). Furt hermore, 6y"f^ + Sy'^^ + Syi^ — 0. As a result, 
Syi{r,r';uj) — Sy[^. Applying the contraction rules (871 to Eq. (E7| and taking into account W — W, we simplify 

sy'id to 

Sy2d - I (^^) ^0 J dr,yo{ri,ri;Cj){strA+tW{r)A_tW{r')striViW{ri))^)^, (ElO) 
which gives the second term of Eq. (90). 
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